{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "f53f2418",
   "metadata": {},
   "source": [
    "# Open Method\n",
    "**강좌**: *수치해석*"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "24713c71",
   "metadata": {},
   "source": [
    "## Newton-Raphson Method\n",
    "\n",
    "### 알고리즘\n",
    "* 함수 $f(x)$ 가 미분 가능하고 연속함수인 경우에 대해서 Taylor Expansion 을 이용하면 다음과 같다.\n",
    "\n",
    "$$\n",
    "f(x_{i+1}) = f(x_{i}) + f'(x_i)(x_{i+1} - x_i) + \\frac{f''(\\xi)}{2!} (x_{i+1} - x_i)^2.\n",
    "$$\n",
    "\n",
    "여기서 $\\xi \\in (x_i, x_{i+1})$ 이다. 선형 근사하고, $x_{i+1}$ 에서 x 축과 교점이 생기는 경우 다음과 같다.\n",
    "\n",
    "$$\n",
    "0 = f(x_{i+1}) = f(x_{i}) + f'(x_i)(x_{i+1} - x_i).\n",
    "$$\n",
    "\n",
    "즉, 아래 식을 반복적으로 해석해서 계산할 수 있다.\n",
    "\n",
    "$$\n",
    "x_{i+1} = x_{i} - \\frac{f(x_i)}{f'(x_i)}\n",
    "$$\n",
    "\n",
    ":::{figure-md} newton-fig\n",
    "<img src=\"https://upload.wikimedia.org/wikipedia/commons/thumb/2/22/Newton_method_scheme.svg/465px-Newton_method_scheme.svg.png\" alt=\"newton-fig\">\n",
    "\n",
    "Newthon Method (from Wikipedia)\n",
    ":::\n",
    "\n",
    "#### 종료 판정 기준\n",
    "Bracketing method와 동일하게 근사 상대 오차 $\\epsilon_a$ 의 크기가 $tol$ 보다 작으면 멈춘다.\n",
    "    \n",
    ":::{note}\n",
    "$\\epsilon_a< tol$ 이면 멈춘다.\n",
    ":::"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a9f3812",
   "metadata": {},
   "source": [
    "### 예제 1\n",
    "$f(x) = e^{-x} - x$ 의 근을 구하시오. 초기 가정은 $x_0 = 0$ 이다.\n",
    "\n",
    "#### By Hand\n",
    "도함수 $f'(x) = -e^{-x} - 1$ 이다.\n",
    "\n",
    "Newton-Raphson 식은 아래와 같다.\n",
    "\n",
    "$$\n",
    "x_{i+1} = x_{i} - \\frac{f(x_i)}{f'(x_i)} = x_i - \\frac{e^{-x} - x}{-e^{-x} - 1}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "1b2a541b",
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "import numpy as np\n",
    "\n",
    "plt.style.use('ggplot')\n",
    "plt.rcParams['figure.dpi'] = 150"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "8f6bab1a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7fe95bf43590>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1wAAAKPCAYAAAB5OMvxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAABcSAAAXEgFnn9JSAABvmElEQVR4nO3deXzU1b3/8feZTPYEwhJIAkJYDIssCii4VFEKVMFq8SpdQO3V+nOpVy/1em2xildvpa1FbS9YK3Vp1bZW0SpYWlHEFdksiyIBAmFJgIQAIXsmc35/TDJJSCaQkMl3MvN6Ph558D3nu+Qz9ZDyzvl+z9dYa60AAAAAAO3O5XQBAAAAABCuCFwAAAAAECQELgAAAAAIEgIXAAAAAAQJgQsAAAAAgoTABQAAAABBQuACAAAAgCAhcAEAAABAkBC4AAAAACBICFwAAAAAECQELgAAAAAIEgIXAAAAAAQJgQsAAAAAgoTABQAAAABB4na6gPaSk5OjTZs2aceOHdq+fbuOHDmi6OhovfTSSy2eV11drb///e/65JNPlJ+fL6/Xq+7du2vo0KGaOXOmunfv3uScbdu2acmSJcrOzpbH41Hfvn01depUTZw4MUifDgAAAEBnFDaB69VXX9W6detadc6xY8f08MMPa8+ePUpJSdHIkSMlSQcOHNDKlSt16aWXNglca9as0YIFC2St1bBhw5ScnKwtW7Zo0aJFys3N1Q033NBunwkAAABA5xY2gSsrK0uZmZkaNGiQBg0apFtuuaXF471er37xi19oz549mjFjhq699lpFRUX59x88eFDx8fGNzikpKdGiRYvk9Xr1ox/9SOPHj5ckHT16VA888ICWLVumsWPHasSIEe3/AQEAAAB0OmETuK6++upWHf/+++9r+/btGj9+vL797W832d+7d+8mfe+9957Kyso0btw4f9iSpJSUFM2aNUuPPfaYli5dSuACAAAAICmCF81YsWKFJOnKK6885XPWr18vSZowYUKTfWPGjFF0dLQ2b96sqqqq9ikSAAAAQKcWNjNcrVFeXq6dO3cqPj5egwcPVnZ2ttatW6eSkhL17NlT48aNU79+/Zqct2fPHknSwIEDm+xzu93q16+fdu7cqby8PGVmZp52nQcOHJC19rSvg9PXs2dPSVJhYaHDlSDUMDbQEsYHAmFsIBDGRugxxigtLa3N50dk4Nq3b5+stUpLS9Nzzz2nf/zjH432/+Uvf9GVV16pWbNm+fvKyspUWloqSc2uXFjXv3PnThUWFrZL4LLWErhCDP89EAhjAy1hfCAQxgYCYWyEj4gMXHXBac+ePdq1a5euvPJKTZ06VXFxcVq7dq2ef/55vfnmm+rVq5emTJkiSaqoqPCfHxsb2+x16/obHtuSOXPmNOmLiYnR/PnzJdX/hgPOc7t9f1VSU1MdrgShhrGBljA+EAhjA4EwNsJPRD7D5fV6JUk1NTW68MILNXv2bPXq1UtdunTRpEmT9L3vfU+S9PrrrztZJgAAAIBOLiJnuOLi4vzbl156aZP9l156qZ577jkdPnxYBw4cUFpaWqNzKisrlZCQ0OS8ysrKJtdvyYIFC1rcX1hYyHRyiKj7LVNBQYHDlSDUMDbQEsYHAmFsIBDGRugxxig9Pb3N50fkDFevXr38281N18bGxqpLly6SfC9HlqSEhAR/yCoqKmr2unX93AoIAAAAQIrQwNWzZ08lJydL8r3M+ERer9f/nFfD2ar+/ftLknJycpqc4/F4tGfPHkVHRysjIyMYZQMAAADoZCIycEnS2LFjJUlffPFFk33Z2dnyeDyKiYlRnz59/P1jxoyRJK1evbrJORs2bFB1dbVGjBihmJiYIFUNAAAAoDOJ2MD1zW9+Uy6XS2+++aZ27drl7z927Jiee+45Sb5nuepWipGkSZMmKT4+XuvWrdNnn33W6JwXX3xRkjR9+vQO+gQAAACRq+71OXzx1dqvjmasE981CDZs2KDXXnvN396+fbuMMRo8eLC/75prrvHPUknS3//+dz333HOKjo5WVlaWYmNjtW3bNpWWlmrAgAGaN2+e4uPjG32f1atX6/HHH5ckDR8+XMnJydq8ebNKS0t1+eWX6/vf/367fab8/HxHBgWa4gFWBMLYQEsYHwiEsdE2Ho9H5eXlqqqqCtt/I7lcvvmQulW10f6MMYqJiVF8fHyjyZWWjj+dRTPCZpXC4uJibd++vVGftbZRX3FxcaP9l19+uTIyMvTWW29px44dqq6uVu/evTVt2jRdeeWVzb5va8KECXrooYe0ZMkSbd++XR6PR3369NHUqVObXfEQAAAAp8/j8ejYsWOKi4tTSkqKP5iEm7oA4PF4HK4kfHm9XlVUVOjYsWPq2rXrKYWu0xE2M1zhiBmu0MFvIhEIYwMtYXwgEMZG6x0/flwul0uJiYlOlxJUBK6OU1paKq/X619MLxCWhQcAAEDYq6qqOuV3nQKnIjY2VlVVVUH/PgQuAAAAhLS6xQ7C9TZCOCMqKqpDFtJg1AIAAABAkBC4AAAAACBICFwAAAAAECQELgAAAAABffDBB+rTp4/Wrl3r7/vnP/+pQYMG6Z133nGwss6BwAUAAAAgoIsvvlgTJ07Uo48+Kkn69NNPdfvtt+uXv/ylJk+e7HB1oS9sXnwMAAAAoGVHjhzR0aNHWzwmOTlZPXv2bNR3//33a8qUKfr1r3+tRYsW6f7779eMGTOCWGn4IHABAAAAEeLZZ5/VggULWjzm2muv1RNPPNGob9iwYbrqqqv085//XPfcc49uvPHG4BUZZrilEAAAAOjkvvrqK91zzz2aMGGCBgwYoNGjR+uuu+5Sfn5+o+P+8z//U7m5uS1+NRfIvvjiC7333nuKjo5Wjx49OupjNXHHHXeoT58+evLJJ5vsW7NmjQYOHKhRo0Zp9+7dHV9cAAQuNGEPH5L3L79XzW8eVs3jDzhdDgAAAFrw/PPPa+rUqfrLX/6iXr16afLkyUpKStKrr76qb37zmyoqKvIf63K55Ha7W/w68QXTu3bt0ve+9z1df/31+uEPf6gFCxaotLS0oz+mJOmee+6R2+3W7373Ox0/ftzfv2PHDn3/+9+Xy+XSH/7wB2VmZjpSX3MIXGiqukp2xd+kTWulrZtkPdVOVwQAAIBm/P3vf9f999+vfv36afny5XrzzTf1u9/9Th988IG+853vKC8vT0899VSbr5+fn6/vfOc7mjp1qu677z7deuutqqmpOa1rno4BAwZo5syZOnr0qBYvXixJKigo0OzZs3X8+HH99re/1dlnn+1IbYEQuNBUz96SqR0a1isVHnS2HgAAADRRWlqq++67T0lJSXrppZd01lln+fdFRUVp7ty5kqQPP/ywTdcvKirSd7/7XY0aNcq/QmFSUpL+4z/+Q08//bQOHnTm34h333234uLi9Mwzz+jAgQO64YYbtGfPHs2fP19f//rXHampJSyagSaMO1rq2UsqOODrOJgvpfV1tigAAIBmWGulcmdubzst8YkyxpzWJV588UUVFhbqjjvuUL9+/Zrs79atm5KSkk66KmEg3bt318qVK5v0/+AHP9APfvCDk57/gx/8QNu2bWvV93zyySd1zjnntHhMRkaGZs+erWeeeUaTJ09WUVGR5syZo+9+97ut+l4dhcCF5vXO8Acue3C/jM51uCAAAIBmlJfKe1do/kO7Ja4nX5YSkk7rGv/4xz8kSQsXLtTChQsDHte/f//T+j5ttW/fPu3cubNV55SXl5/ScbfccosWL16soqIiXXfddfrRj37UlhI7BIELzTK9MmS1wdc4lOdsMQAAAGjiyy+/VFxcnK688soWj2t4q2FH+vvf/x6U61pr9dBDD/lmNyW53aEdaUK7Ojind4Z/0x7Kb+FAAAAAdLSqqiodP35c/fr1a/LOrHA3b948LV26VJMnT9bnn3+uV155RbfddpsGDhzodGnNInChWb4ZrloH9ztZCgAAQGDxib7b8zqb+MTTOj0mJkZxcXHKy8tTWVmZEhIS2qmw0Pb0009r8eLFOuecc/TUU0/pD3/4g/7nf/5Hjz32mBYtWuR0ec1ilUI0r8EMl4oKZasqnasFAAAgAGOMTEJS5/s6zQUzJOniiy+Wx+PRfffd1+TZJ2utPvzwQ61du/a0v0+oePPNN/Xwww8rMzNTL7zwguLj43X99derV69eevPNN/Xll186XWKzCFxoXo9UKarBBGjdioUAAAAICXPnzlVKSopee+01jR8/XrNnz9att96qq6++Wuecc46+/e1vO/aC4vb26aef6u6771b37t310ksvqUePHpKk+Ph43XHHHbLW6uc//7nDVTaPwIVmGVeUlJpW38FthQAAACFl8ODB+uc//6lZs2YpMTFRH330kd577z0dOnRIY8eO1S9/+UuNHz/e6TJPW3Z2tm666Sa5XC698MILyszMbLR/1qxZSktL04oVK7Ru3TpnimwBz3AhsN4Z0oF9kiR7MF+nP/ENAACA9tSnT5+QndlpL1lZWS3eLhgXF6f169d3YEWtwwwXAjK90usbLA0PAAAAtBqBC4H17uPftNxSCAAAALQagQsBNZ7h4l1cAAAAQGsRuBBYw6Xhjx2RrShzrhYAAACgEyJwIbCUHlJMTH37ILNcAAAAQGsQuBCQcbmk1PrbCi0LZwAAAACtQuBCyxreVniQwAUAAAC0BoELLTIELgAAAKDNCFxoWa/6wMUthQAAAEDrELjQItMgcPHyYwAA4ARjjCTJWutwJQgndeOpbnwFC4ELLUtrELhKjsuWHneuFgAAELHcbreqq6udLgNhpLq6Wm63O+jfh8CFliWnSHHx9W2e4wIAAA6IiYlRRUUFs1xoF9ZaVVRUKKbhK5CCJPiRDp2aMcb3HNeenZIkezBPZuAQh6sCAACRJj4+XpWVlSopKVF8fLyioqKcLiko6gIlwTJ4ampqVF5eLq/Xq/j4+JOfcJoIXDgp0ztDtjZw6eB+Z4sBAAARyRijrl27qry8XMeOHQvbQOJy+W5A83q9DlcSvowxiouLU9euXYP+/JZE4MKpSOtTv32AwAUAAJzhcrmUmJioxMTEsA1cqampkqSCggKHKwlfHRGyGiJw4eR61wcuywwXAAAIAR39j+aOUve5wvXzRSIWzcBJmbS+9Y2DebJMcQMAAACnhMCFk+vdYGn46iqpiCluAAAA4FQQuHBSJi5e6tazvoPnuAAAAIBTQuDCqUnjOS4AAACgtQhcOCWmd8OVCvc5VwgAAADQiRC4cGoaLJxhuaUQAAAAOCUELpwSw7u4AAAAgFYjcOHUNAxcRw/LVpQ5VwsAAADQSRC4cGq69ZRiYurbB/OcqwUAAADoJNxOF9BecnJytGnTJu3YsUPbt2/XkSNHFB0drZdeeumUr/Hwww9r8+bNkqTf/e53SklJafa4bdu2acmSJcrOzpbH41Hfvn01depUTZw4sR0+SWgyLpfUq4+0b5ck33Ncpv9gh6sCAAAAQlvYBK5XX31V69ata/P577//vjZv3ixjjKy1AY9bs2aNFixYIGuthg0bpuTkZG3ZskWLFi1Sbm6ubrjhhjbXEOpMWh/Z2sDFc1wAAADAyYVN4MrKylJmZqYGDRqkQYMG6ZZbbjnlc4uLi/XHP/5Ro0ePVl5engoKCpo9rqSkRIsWLZLX69WPfvQjjR8/XpJ09OhRPfDAA1q2bJnGjh2rESNGtMtnCjlpLA0PAAAAtEbYPMN19dVX67rrrtPYsWMD3goYyPPPP6+KigrddNNNLR733nvvqaysTOPGjfOHLUlKSUnRrFmzJElLly5tde2dRoN3cbE0PAAAAHByYRO42upf//qXPvroI82YMUNpaWktHrt+/XpJ0oQJE5rsGzNmjKKjo7V582ZVVVUFpVanNVoa/tB+Wa/XuWIAAACATiCiA1dlZaWeeeYZ9enTR1ddddVJj9+zZ48kaeDAgU32ud1u9evXT9XV1crLC9MV/BoGrqoq6chh52oBAAAAOoGIDlx/+ctfVFBQoJtvvllud8uPs5WVlam0tFSS1L1792aPqesvLCxs30JDhIlLkFIafHae4wIAAABaFDaLZrRWTk6O3n77bV1yySU666yzTnp8RUWFfzs2NrbZY+r6Gx7bkjlz5jTpi4mJ0fz58yVJPXv2PKXrdKSiMwao6miRJCmx5JgSU1Mdrqhj1AXy1Aj5vDh1jA20hPGBQBgbCISxEX4icobL6/Xq6aefVmJiombPnu10OZ1KVJ9+/u2avFwHKwEAAABCX0TOcC1btky7du3Srbfeqi5dupzSOXFxcf7tyspKJSQkNDmmsrKyybEtWbBgQYv7CwsLW3wnmBO8XXv4t8t27VBlgCX0w03db5kCvTIAkYuxgZYwPhAIYwOBMDZCjzFG6enpbT4/IgPX+vXrZYzRqlWr9MEHHzTad/ToUUnSY489JrfbrW9/+9saOnSoEhISlJCQoLKyMhUVFTUbuIqKfLfaheKtgO3FpPWRPwLm8wwXAAAA0JKIDFySZK3V1q1bA+7Pzs6W5Hspcp3+/ftr69atysnJUd++fRsd7/F4tGfPHkVHRysjIyM4RYeC9DPqt48eli0vk4lvGj4BAAAARGjgmjdvXsB9d9xxhwoKCvS73/2uyQuUx4wZo61bt2r16tW6+OKLG+3bsGGDqqurdc455ygmJiYIVYeIbj2l2DipsnZhkAP7pAFZztYEAAAAhKiIXDSjrSZNmqT4+HitW7dOn332mb//2LFjevHFFyVJ06dPd6q8DmFcLimtfnbP5u11sBoAAAAgtIXNDNeGDRv02muvNerzeDyaO3euv33NNddozJgxbf4eSUlJuu222/T4449rwYIFGj58uJKTk7V582aVlpbq8ssv18iRI9t8/c7CpPeVzd3ha+QTuAAAAIBAwiZwFRcXa/v27Y36rLWN+ho+j9VWEyZM0EMPPaQlS5Zo+/bt8ng86tOnj6ZOnapLL730tK/fKTSc4eLlxwAAAEBAYRO4Jk6cqIkTJ572dRYuXHjSY4YOHaqf/OQnp/29OiuT0a/BSoXMcAEAAACB8AwXWi+9wQqNBQdlq6ucqwUAAAAIYQQutF5quhRVOzlqvdLB/c7WAwAAAIQoAhdazURFSb3q37ZteQEyAAAA0CwCF9omo8ELkHmOCwAAAGgWgQttYtIbBi5muAAAAIDmELjQNg2XhmeGCwAAAGgWgQttYjL61TcO7petqXGuGAAAACBEEbjQNr0zJGN82x6PVHjQ2XoAAACAEETgQpuYmFipZ+/6Dm4rBAAAAJogcKHtGj3HxcIZAAAAwIkIXGgz02hp+D3OFQIAAACEKAIX2q7B0vDMcAEAAABNEbjQZqbBLYU6sE/WWueKAQAAAEIQgQtt1/DlxxXl0pFC52oBAAAAQhCBC21mEhKllO71HXmsVAgAAAA0RODC6WnwAmSbx8IZAAAAQEMELpwWk9G/vkHgAgAAABohcOH0NFganhkuAAAAoDECF06LaXBLofL2slIhAAAA0ACBC6enYeCqLJeKCpyrBQAAAAgxBC6cFhOfIHVPre/Yn+tcMQAAAECIIXDh9LFSIQAAANAsAhdOW+PnuAhcAAAAQB0CF05fnwYzXPsJXAAAAEAdAhdOW6MZrgN7Zb1e54oBAAAAQgiBC6cvvf5dXKqqkgoPOlcLAAAAEEIIXDhtJjZOSk2r78hjpUIAAABAInChvWTwHBcAAABwIgIX2oXJaHBbYd5e5woBAAAAQgiBC+0jo79/03JLIQAAACCJwIV20nilwn2yNTXOFQMAAACECAIX2kd6X8nUDiePRyrId7YeAAAAIAQQuNAuTHSM1Cu9viOPhTMAAAAAAhfaT4OFM1ipEAAAACBwoR2ZPvULZ2g/C2cAAAAABC60m4aBy+7f7VwhAAAAQIggcKH99M2s3z6YL1tV6VgpAAAAQCggcKH99EqXomN829Yr5fMCZAAAAEQ2AhfajXFFSQ3ex2X37XauGAAAACAEELjQrhotnLGPhTMAAAAQ2QhcaF8NnuNi4QwAAABEOgIX2lXjGa7djtUBAAAAhAICF9pXw5UKjx+TLT7iWCkAAACA0whcaFemS4qU3LW+g+e4AAAAEMEIXGh/DZ/j4rZCAAAARDACF9qd6ZNZ39jPDBcAAAAil9vpAtpLTk6ONm3apB07dmj79u06cuSIoqOj9dJLLzU51uv1atu2bVq/fr2+/PJLHTp0SGVlZerRo4dGjhypq6++Wr169Qr4vbZt26YlS5YoOztbHo9Hffv21dSpUzVx4sQgfsJOhBkuAAAAQFIYBa5XX31V69atO6VjDx06pAcffFCS1L17d2VlZcnlcmnHjh1asWKFPv74Y/34xz/W0KFDm5y7Zs0aLViwQNZaDRs2TMnJydqyZYsWLVqk3Nxc3XDDDe36uToj07e/bF0jf69sTY1MVJSTJQEAAACOCJvAlZWVpczMTA0aNEiDBg3SLbfc0uLxo0eP1re+9S0NHz7c31ddXa1nnnlG77//vn7961/r17/+tdzu+v+JSkpKtGjRInm9Xv3oRz/S+PHjJUlHjx7VAw88oGXLlmns2LEaMWJEcD5kZ5F+hmRckvVK1VXSoXwpva/TVQEAAAAdLmye4br66qt13XXXaezYsUpJSWnx2LS0NM2dO7dR2JKk6Oho3XzzzUpISFBhYaGys7Mb7X/vvfdUVlamcePG+cOWJKWkpGjWrFmSpKVLl7bPB+rETEys1Du9voMXIAMAACBChU3gai8xMTFKT/eFhaKiokb71q9fL0maMGFCk/PGjBmj6Ohobd68WVVVVcEvNMQ1XDiD57gAAAAQqQhcJ/B6vSosLJSkJjNle/bskSQNHDiwyXlut1v9+vVTdXW18vLygl5nyOvb379pWakQAAAAEYrAdYKPP/5Yx44dU5cuXTRkyBB/f1lZmUpLSyX5FtpoTl1/XWCLZKbBSoVihgsAAAARKmwWzWgPhYWFev755yVJ1113naKjo/37Kioq/NuxsbHNnl/X3/DYlsyZM6dJX0xMjObPny9J6tmz5yldJxR5Ro2VP3YWHlSPxAS5EhKdLOm01C2ekpqa6nAlCDWMDbSE8YFAGBsIhLERfpjhqlVRUaHHHntMx48f17nnnqspU6Y4XVKnFtUrXaZBwPLs3uFgNQAAAIAzmOGS5PF49Ktf/Uo5OTkaOnSo7rrrribHxMXF+bcrKyuVkJDQ5JjKysomx7ZkwYIFLe4vLCyUtbbFY0KZ7dNf2v6lJOnI5g1ypWY4XFHb1f2WqaCgwOFKEGoYG2gJ4wOBMDYQCGMj9Bhj/IvqtUXEz3B5vV795je/0caNG9W/f3/993//t2JiYpocl5CQ4A9ZJ65eWKeuvzPfCtiezBkNFhfZu8u5QgAAAACHRHzgWrx4sT799FOlp6fr/vvvV2Ji4OeM+vf3rbyXk5PTZJ/H49GePXsUHR2tjIzOO5PTrs4Y4N+0BC4AAABEoIgOXC+//LJWrFihnj176qc//am6du3a4vFjxoyRJK1evbrJvg0bNqi6ulojRoxodoYsEjWa4dqfK+vxOFcMAAAA4ICIDVxLly7VG2+8oZSUFP30pz89pdsAJ02apPj4eK1bt06fffaZv//YsWN68cUXJUnTp08PWs2dTsYZUlSUb9tTLR3c72w9AAAAQAcLm0UzNmzYoNdee61Rn8fj0dy5c/3ta665RmPGjNHu3bv1xz/+UZLUq1cvLVmypNlrTpo0SUOHDvW3k5KSdNttt+nxxx/XggULNHz4cCUnJ2vz5s0qLS3V5ZdfrpEjRwbh03VOJjpGSusr1b742O7NkenT/yRnAQAAAOEjbAJXcXGxtm/f3qjPWtuor7i4WJJUWlrqX/0vOztb2dnZzV7zrLPOahS4JGnChAl66KGHtGTJEm3fvl0ej0d9+vTR1KlTdemll7bnRwoL5oyBsrWBS3t3SxMcLQcAAADoUMZ25nXHw1x+fn6nXhZekrz/fEP2r8/6GsNGK2rOw84W1EYs0YpAGBtoCeMDgTA2EAhjI/SwLDxCmmmwUqH27ur0ARIAAABoDQIXgqtvg8BVUiwdbf4dZgAAAEA4InAhqExyF6lbgxUg9zZ9hxkAAAAQrghcCD5egAwAAIAIReBC0J34HBcAAAAQKQhcCDpzxkD/NjNcAAAAiCQELgRfwxmugnzZijLnagEAAAA6EIELwdeztxQX79u2Vtq329FyAAAAgI5C4ELQGZer0fLw3FYIAACASEHgQocw/QfVN3J3OFcIAAAA0IEIXOgY/RosnJHLu7gAAAAQGQhc6BCmX4MZrvw9stVVzhUDAAAAdBACFzpG+hlSdIxvu6ZG2p/rbD0AAABAByBwoUOYqCipb6a/bXN3OlcMAAAA0EEIXOgwjRbO2EPgAgAAQPgjcKHjNHiOixkuAAAARAICFzpMo4Uz9u+W9XicKwYAAADoAAQudJw+/aQot2/b45Hy9jhbDwAAABBkBC50GOOOlvr097ctz3EBAAAgzBG40KFYOAMAAACRhMCFjtVvoH+ThTMAAAAQ7ghc6FCm/+D6xr5dsjU1zhUDAAAABBmBCx2rT3/JVTvsqqqkA/udrQcAAAAIIgIXOpSJiZXSz/C3WTgDAAAA4YzAhQ7X6H1cuTucKwQAAAAIMgIXOl6D57hYOAMAAADhjMCFDmcyGyycsWcnC2cAAAAgbBG40PH6DmiwcEallL/X2XoAAACAICFwocOZ2Fgpo7+/bXdvd7AaAAAAIHgIXHCEGXBmfYOFMwAAABCmCFxwRoPnuOwuZrgAAAAQnghccITJbDDDtW+3bHW1c8UAAAAAQULggjMy+kvRMb7tGo+0b7ej5QAAAADBQOCCI4zbLZ0xwN9m4QwAAACEIwIXHGMavABZBC4AAACEIQIXnNPgOS7LSoUAAAAIQwQuOKbR0vB5e2UrK5wrBgAAAAgCAhec07uPFBfv27ZeKXens/UAAAAA7YzABccYl0tq8BwXC2cAAAAg3BC44CiTycIZAAAACF8ELjiq4QuQmeECAABAuCFwwVkNl4YvOCBbety5WgAAAIB2RuCCs3r2lpK71rd3ZTtXCwAAANDOCFxwlDFGGpDlb9scAhcAAADCB4ELjjMNA9eubQ5WAgAAALQvAhccZwYOqW/kZMta61wxAAAAQDsicMF5mWdKxvi2y0qkg3nO1gMAAAC0E7fTBbSXnJwcbdq0STt27ND27dt15MgRRUdH66WXXmrxvFWrVmn58uXat2+f3G63srKyNGPGDA0ZMiTgOdu2bdOSJUuUnZ0tj8ejvn37aurUqZo4cWI7f6rIYBISpbS+Uv5eSZLdlS2T1sfhqgAAAIDTFzaB69VXX9W6detadc4LL7ygZcuWKSYmRqNGjVJ1dbU2bdqkjRs3as6cOTrvvPOanLNmzRotWLBA1loNGzZMycnJ2rJlixYtWqTc3FzdcMMN7fWRIooZmCVbG7iUs006/1JnCwIAAADaQdgErqysLGVmZmrQoEEaNGiQbrnllhaP37Jli5YtW6bk5GQ98sgjSk9PlyRlZ2dr3rx5WrRokYYPH66kpCT/OSUlJVq0aJG8Xq9+9KMfafz48ZKko0eP6oEHHtCyZcs0duxYjRgxIngfNFwNHCJ9/K4k3wwXAAAAEA7C5hmuq6++Wtddd53Gjh2rlJSUkx7/1ltvSZJmzJjhD1uSL7hNnjxZZWVlWrlyZaNz3nvvPZWVlWncuHH+sCVJKSkpmjVrliRp6dKl7fBpIo8Z0OAWzn27ZKsqnSsGAAAAaCdhE7hao6qqSlu2bJEkTZgwocn+ur7169c36q9rN3fOmDFjFB0drc2bN6uqqqq9Sw5/Gf2kmFjfdk2NtGens/UAAAAA7SAiA1deXp6qq6vVpUsX9ejRo8n+AQMGSJJyc3Mb9e/Zs0eSNHDgwCbnuN1u9evXT9XV1crLY5W91jJRUb7VCmvZHN7HBQAAgM4vIgNXYWGhJDUbtiQpLi5OiYmJKi0tVXl5uSSprKxMpaWlkqTu3bs3e15df9310ToNX4CsHJ7jAgAAQOcXNotmtEZFRYUkKSYmJuAxsbGxKi0tVUVFheLj4/3n1O0LdE7D65/MnDlzmvTFxMRo/vz5kqSePXue0nXCRcU55+roP5ZIkkzuDqWmpjpcUT232/dXJZRqQmhgbKAljA8EwthAIIyN8BORM1zWWkmSqXvZbgvHoONEZ53l3/YWHlRNUYGD1QAAAACnLyJnuOLj4yVJlZWBV8KrW/giLi6u0Z915yUkJDQ5p+56DY9tyYIFC1rcX1hYGHnBr3tPqch3S+bhtZ/KjDnf4YJ86n7LVFBACERjjA20hPGBQBgbCISxEXqMMY1WNW+tiJzhqrtV7/Dhw83ur6ioUGlpqRITE/3hLCEhwR+yioqKmj2vrj/SbgVsVw2e47I5XzlYCAAAAHD6IjJwZWRkKDo6WsXFxc2Grl27dkmS+vXr16i/f//+kqScnJwm53g8Hu3Zs0fR0dHKyMgIQtWRwQwa5t+2OwlcAAAA6NwiMnDFxMRoxIgRkqTVq1c32V/XN3bs2Eb9Y8aMCXjOhg0bVF1drREjRrS4GAdaZgYNrW/s3iFbXe1cMQAAAMBpisjAJUnTpk2TJC1ZskT5+fn+/uzsbK1YsULx8fG67LLLGp0zadIkxcfHa926dfrss8/8/ceOHdOLL74oSZo+fXoHVB/G+g2UomsDq6eaFyADAACgUwubRTM2bNig1157rVGfx+PR3Llz/e1rrrnGP0s1atQoXXHFFXr77bd17733auTIkaqpqdGmTZvk9Xp15513KikpqdH1kpKSdNttt+nxxx/XggULNHz4cCUnJ2vz5s0qLS3V5ZdfrpEjRwb/w4Yx446WMgdL27+UJNmdWxvPegEAAACdSNgEruLiYm3fvr1Rn7W2UV9xcXGj/TfeeKMyMzO1fPlybd68WVFRURoxYoSuueYaDR3a/D/yJ0yYoIceekhLlizR9u3b5fF41KdPH02dOlWXXnpp+3+wCGQGDZP1By6e4wIAAEDnZWzErTveeeTn50fesvCS7MY18v7fI75GlxS5HnuhxXemdQSWaEUgjA20hPGBQBgbCISxEXpYFh7hZ2CD2cXio1LhQcdKAQAAAE4HgQshxyR3kXr38bftzq0OVgMAAAC0HYELIckMbjDLtYPABQAAgM6JwIXQxAuQAQAAEAYIXAhJjZaC358rW17mXDEAAABAGxG4EJrS+koJte9Bs1bK2eZsPQAAAEAbELgQkozLJTWY5WLhDAAAAHRGBC6ErIa3FVoWzgAAAEAnROBCyDKD6xfOUE62bE2Nc8UAAAAAbUDgQujKzJKi3L7tynJpb46z9QAAAACtROBCyDKxsVLmYH/bbv/SwWoAAACA1iNwIaSZwcP923b7Fw5WAgAAALQegQshzZx5Vn1j+5ey1jpXDAAAANBKBC6EtsHDJGN82yXF0oF9ztYDAAAAtAKBCyHNJCZJffr729xWCAAAgM6EwIWQZ86sf45L2QQuAAAAdB4ELoS+M0f4N1mpEAAAAJ0JgQshr9EMV1GB7OFDzhUDAAAAtAKBCyHPpHSXUtP8bZ7jAgAAQGdB4EKnYLIaLw8PAAAAdAYELnQODd7HxXNcAAAA6CwIXOgUGj3Hlb9X9vgx54oBAAAAThGBC51DarrUtXt9m+XhAQAA0AkQuNApGGMaPcdls7c4WA0AAABwaghc6DyGjvRv2m2bHSwEAAAAODUELnQaJqs+cGl/Ls9xAQAAIOQRuNB59M7gOS4AAAB0KgQudBrGGJkhDW8r3ORgNQAAAMDJEbjQuQwZ4d+021g4AwAAAKGNwIVOxTRYOEN5e2SLjzpWCwAAAHAyBC50LqnpUkqP+jbLwwMAACCEEbjQqfie4+K2QgAAAHQOBC50PkN4HxcAAAA6BwIXOp2GKxUqf69s8RHnigEAAABaQOBC55OaJnXv6W/abbyPCwAAAKGJwIVOxxgjk9Vglov3cQEAACBEEbjQOTVYHt5uJXABAAAgNBG40CmZoaPrG4fyZA8fcq4YAAAAIAACFzol0yNV6t3H37ZbNzpYDQAAANA8Ahc6LTNsVH2D2woBAAAQgghc6LTMsPrbCu3Wf8la62A1AAAAQFMELnReQ0ZJxvi2jx+T9uc6Ww8AAABwAgIXOi2TmCT1G+Rv2694jgsAAAChhcCFTs0Mb3Bb4ZcELgAAAIQWAhc6tUbLw2dvkfV4nCsGAAAAOAGBC53b4GGSO9q3XVkh7cp2th4AAACgAQIXOjUTEyudOdzf5n1cAAAACCVupwtwWnZ2tt58801t27ZNJSUliouL04ABAzRlyhRNmDCh2XNWrVql5cuXa9++fXK73crKytKMGTM0ZMiQDq4ekmSGjvIHLbt1o/TN7zhcEQAAAOAT0YHr008/1RNPPCFrrQYNGqSzzjpLR44c0RdffKEtW7boqquu0ve+971G57zwwgtatmyZYmJiNGrUKFVXV2vTpk3auHGj5syZo/POO8+hTxO5zLCzZV//o6+xa5tseZlMfIKzRQEAAACK4MBVU1Oj3//+97LW6u6779YFF1zg35edna2HHnpIb775piZNmqS0tDRJ0pYtW7Rs2TIlJyfrkUceUXp6uv/4efPmadGiRRo+fLiSkpIc+UwRq/9AKTFZKj0u1dRI2zZJZzc/OwkAAAB0pIh9hmv//v0qLi5Wnz59GoUtScrKytLo0aNlrVVOTo6//6233pIkzZgxwx+26o6fPHmyysrKtHLlyo75APAzriiZYQ2Wh//icwerAQAAAOpFbOCKjo4+pePqZquqqqq0ZcsWSWr22a66vvXr17dThWiVs87xbxK4AAAAECoiNnD17t1bvXv31v79+/XJJ5802pedna2NGzeqV69eGj7ctwJeXl6eqqur1aVLF/Xo0aPJ9QYMGCBJys3NDX7xaMIMrw9cKjggeyjPuWIAAACAWhH7DJfL5dLtt9+un//853riiSf01ltvqXfv3jpy5Ii++uorDR48WHfeeafcbt//RIWFhZLUbNiSpLi4OCUmJqq0tFTl5eWKj48/aQ1z5sxp0hcTE6P58+dLknr27NnWjxd5UlNV2G+gPHt8t4Am5m5X4lmjT3LSqasbB6mpqe12TYQHxgZawvhAIIwNBMLYCD8RO8MlScOGDdO8efPUq1cv7dy5U5988om2bt2quLg4jRw5Ut26dfMfW1FRIckXiAKJjY1tdCw6VszZ4/3bVZ9/5mAlAAAAgE/EznBJ0kcffaSnnnpKZ555pu6++2717dtXR44c0VtvvaUlS5Zoy5Ytmjdvntxut6y1kiRjTMDr1R1zqhYsWNDi/sLCwlZfM5LZgUP925Wb1utQfp6M+9Se1TuZut8yFRQUtMv1ED4YG2gJ4wOBMDYQCGMj9BhjGi2Y11oRO8OVn5+vhQsXqkuXLrrvvvs0ePBgxcXFKT09XbfccovGjh2r7Oxsvf/++5Lkv0WwsrIy4DWrqqok+W4vhAPOHC5F185AVpZLO79yth4AAABEvIgNXB9//LFqamo0evToZgPS+eefL0n64osvJNU/T3X48OFmr1dRUaHS0lIlJiae0vNbaH8mJlbKOsvftl9scLAaAAAAIIIDV1FRkSQpISGh2f11/SUlJZKkjIwMRUdHq7i4uNnQtWvXLklSv379glEuTpE5a4x/m+XhAQAA4LSIDVwpKSmSpJ07dza7f8eOHZLq76ONiYnRiBEjJEmrV69ucnxd39ixY9u7VLSCafA+Lu3JkS0+4lwxAAAAiHgRG7jGjRsnSdq6dav++c9/NtqXnZ2tZcuWSWr8kuNp06ZJkpYsWaL8/PxGx69YsULx8fG67LLLgl06WpJ+htStfjl9+8W/nKsFAAAAES9iVykcOHCgrrzySr311ltavHix/vGPf6hPnz46cuSIsrOzZa3V17/+dY0aNcp/zqhRo3TFFVfo7bff1r333quRI0eqpqZGmzZtktfr1Z133qmkpCQHPxWMMTIjxsh+WBuiN6+Tzr/U2aIAAAAQsSI2cEnS7NmzNWTIEL3zzjvKyclRXl6e4uLiNHz4cE2aNEkXXXRRk3NuvPFGZWZmavny5dq8ebOioqI0YsQIXXPNNRo6dGgz3wUdzYwc5w9c9osNsjU1MlFRDlcFAACASBTRgUuSzjvvPJ133nmtOmfixImaOHFicArC6Rs2WnK7JY9HKiv1LQ/fYPVCAAAAoKNE7DNcCF8mLl7KGuFv2y3rHKwGAAAAkYzAhbBkRtavFmk3EbgAAADgDAIXwpIZeW59Y3+u7OEC54oBAABAxGr3wHXnnXfqjTfe0LFjx9r70sApM70zpF7p/rbdzCwXAAAAOl67B65Dhw7pT3/6k2677TYtWLBAmzZtau9vAZwSM3Kcf9tuWe9gJQAAAIhU7R64vvWtb6l79+6qqanRZ599pv/93//VnXfeqb/97W/MeqFDmVH1gUtbN8pWVzlXDAAAACKSsdba9r6o1+vV559/rnfeeUcbN26U1+uVJEVFRencc8/VpEmTGr1QGM3Lz89XEP7zRAxbXS3vf35PqqyQJLnuelBmxNiTnNW81NRUSVJBAc+CoTHGBlrC+EAgjA0EwtgIPcYYpaenn/zAAILyHi6Xy6WxY8dq7NixKioq0nvvvaeVK1eqsLBQq1ev1urVq9WrVy99/etf18SJE9W1a9dglIEIZ6Kjfe/k+tdnknyrFbY1cAEAAABtEZQZruZYa7Vx40atWLFCGzZsUE1NjSRmvVrCDNfp837wD9k/LvQ1evSS69FnZIxp9XX4bRMCYWygJYwPBMLYQCCMjdATkjNczTHG6Oyzz9bZZ5+toqIi/frXv9bWrVtVU1Pjn/VKS0vTlVdeqcsuu0wuFyvW4/SZUePkj6yHD0n7d0t9BzhYEQAAACJJh6aawsJCvfLKK5o7d662bt3q78/MzJTL5dKBAwf0zDPPaO7cuSouLu7I0hCmTEoPKfNMf9v+a42D1QAAACDSBH2Gy+v1at26dXr33Xe1adMm/wIaSUlJmjhxoiZPnqy0tDQdPXpU//znP7Vs2TLl5OTo5Zdf1q233hrs8hABzNnjZXdvlyTZjWuk6TMdrggAAACRImiB69ChQ3r33Xf1/vvv6+jRo/7+rKwsTZkyReeff77c7vpvn5KSouuuu05jx47VT37yE33++efBKg0Rxow+T/aNF32N3dtljx72zXwBAAAAQdbugWv16tVasWKFtmzZ4l/wIT4+XhdddJGmTJmifv36tXj+oEGDlJKS0iikAaelT3+pZ2+p8KAkyW5cK3PJNxwuCgAAAJGg3QPX448/7t/OzMzU5MmTddFFFykuLu7Ui3J32FoeiADGGN8s17tvSZLsvz6TCFwAAADoAO2ebKKjo3XBBRdoypQpGjx4cJuusXDhwnauCpGuYeDSVxtlK8pl4uKdLQoAAABhr90D19NPP63ExMT2vixwes48S0pIlMpKJY9H+uJzaewFTlcFAACAMNfuy8ITthCKjNstM2Kcv203fuZgNQAAAIgUvF0YkePs8f5Nu2mdbE2Ng8UAAAAgEhC4EDHMiDFSVO1dtKXHpR1fOlsQAAAAwh6BCxHDxCdIw0b523bDpw5WAwAAgEhA4EJEMeec79+2Gz6V9XodrAYAAADhjsCFiGLOHi+Z2mF/9LC0e7uzBQEAACCsEbgQUUyXFOnM4f42txUCAAAgmAhciDhmTMPbCj+RtdbBagAAABDOCFyIOOacCfWNggPS/t2O1QIAAIDwRuBCxDHdU6UBWf42txUCAAAgWAhciEiNbyskcAEAACA4CFyISA0Dl/bnyh7Y71wxAAAACFsELkQk0ytD6pvpb9sNnzhXDAAAAMIWgQsRy4y5wL9t13/sYCUAAAAIVwQuRCwz7sL6xp4c2YN5zhUDAACAsETgQsQy6WdIffr723bdRw5WAwAAgHBE4EJEM+d+zb9t137oYCUAAAAIRwQuRDQz7qL6xv5c2fy9zhUDAACAsEPgQkQzvTOkfgP9bbuW2woBAADQfghciHhmXIPbCnmOCwAAAO2IwIWI12i1wvy9svtznSsGAAAAYYXAhYhnUtOkzDP9bWa5AAAA0F4IXIAaL55h13woa62D1QAAACBcELgAnbBa4aE8KXeHc8UAAAAgbBC4AEmmR6p05nB/2362ysFqAAAAEC4IXEAtM36if9uu/VDWW+NcMQAAAAgLBC6glhl7gRQV5WscOyJ9tcnZggAAANDpEbiAWiapizRirL9tV3NbIQAAAE4PgQtowIy/xL9tP/9UtqrSwWoAAADQ2bmdLiAUHD16VG+88YY+//xzFRYWKiYmRr169dLIkSM1a9asJsevWrVKy5cv1759++R2u5WVlaUZM2ZoyJAhDlSP9mRGnScbGy9VlksV5bIb18qce9HJTwQAAACaEfEzXNnZ2frP//xPvf3224qKitK4ceOUlZWlkpISLV26tMnxL7zwghYuXKi9e/dq5MiRGjx4sDZt2qQHH3xQa9asceAToD2Z2FiZMRP8bfvZ+84VAwAAgE4vome4ioqK9Oijj6q6ulr33HOPzjvvvEb7d+xo/C6mLVu2aNmyZUpOTtYjjzyi9PR0Sb7QNm/ePC1atEjDhw9XUlJSh30GtD8zfqLspyt9jS0bZEuKpdRUZ4sCAABApxTRM1wvv/yySktLNWvWrCZhS5IGDx7cqP3WW29JkmbMmOEPW5KUlZWlyZMnq6ysTCtXrgxu0Qi+oaOkLim+7RqP7LqPHC0HAAAAnVfEBq6SkhJ9+umnSkhI0GWXXXbS46uqqrRlyxZJ0oQJE5rsr+tbv359+xaKDmeiomTOa7B4xifvOVgNAAAAOrOIvaVw27Ztqq6u1siRI+V2u7V69Wp99dVX8ng86tOnj84//3ylpKT4j8/Ly1N1dbW6dOmiHj16NLnegAEDJEm5ubkd9REQROaCy2RX/M3X2JUtz97dcp+R6WhNAAAA6HwiNnDt3btXktS1a1c98MADys7ObrT/5Zdf1u23367zzz9fklRYWChJzYYtSYqLi1NiYqJKS0tVXl6u+Pj4IFaPYDNnDJDOGCDt3SVJKl/5tpKvv93hqgAAANDZRGzgKi0tlSR98MEHcrvduvXWWzVu3DhVVFRo+fLlWrp0qX7zm98oIyND/fv3V0VFhSQpJiYm4DVjY2NVWlqqioqKUwpcc+bMadIXExOj+fPnS5J69uzZlo+GdlI6+Zs6/uyTkqTy95cr5YY7lMriGTiB2+37McrYQHMYHwiEsYFAGBvhJ2Kf4fJ6vZKkmpoa3XDDDbrsssvUpUsX9erVS9dff70mTJggj8ejv/3Nd1uZtVaSZIwJeM26YxAe4i6eIkVFSZK8RQWq/BfL/gMAAKB1InaGq24GyhijSy65pMn+Sy+9VKtXr9aXX37Z6PjKysqA16yqqpLku73wVCxYsKDF/YWFhYQ4p40YK230Ba2SFW+pOiPT2XoQcup+A1lQUOBwJQhFjA8EwthAIIyN0GOMabRCeWtF7AxX3WBOSUlRdHR0wP3Hjh2TVH973+HDh5u9XkVFhUpLS5WYmMjzW2HEdcEk/3bFZ6tky0odrAYAAACdTcQGrrpVBUtLS5udRSopKZFUP1uVkZGh6OhoFRcXNxu6du3yLa7Qr1+/YJUMJ4waJyUl+7arqmTXfehsPQAAAOhUIjZw9evXT7169VJVVZW2b9/eZP8XX3whSRo4cKAk32IWI0aMkCStXr26yfF1fWPHjg1WyXCAcUc3fifXx+86WA0AAAA6m4gNXJJ01VVXSZKee+45FRcX+/tzcnK0dOlSSdLkyZP9/dOmTZMkLVmyRPn5+f7+7OxsrVixQvHx8af0EmV0LubC+tsKlbNNdj/vWgMAAMCpidhFMyRp0qRJ2rx5s1avXq27775bWVlZqqys1LZt2+TxeDRp0iRNmDDBf/yoUaN0xRVX6O2339a9996rkSNHqqamRps2bZLX69Wdd96ppKQkBz8RgsH0GyT3oKHy7PxKkmQ//KfMt3/gcFUAAADoDCI6cLlcLt19991655139N577/lvIxw0aJAmT56siy++uMk5N954ozIzM7V8+XJt3rxZUVFRGjFihK655hoNHTq0oz8COkjC5CtVXBe4Vr8ve80NMtGB38kGAAAASJKxrDsesvLz81kWPkT0SIhXwb9fKVvpewG2+cE9cp3XNJAj8rB8L1rC+EAgjA0EwtgIPSwLD3QAV2KS4ho8y2U/esfBagAAANBZELiAUxQ/+cr6xtaNsofyAx8MAAAAiMAFnLLooaOk9DP8bfvxCgerAQAAQGdA4AJOkTFG5qL61wTYj9+V9XgcrAgAAAChjsAFtII5/1IpqnZxz2NF0qY1zhYEAACAkEbgAlrBJHeVGXuBv+19/+8OVgMAAIBQR+ACWslccnl9Y+tG2QP7nSsGAAAAIY3ABbTWmcOlPv39TbtquYPFAAAAIJQRuIBWMsY0muWyn7wrW1npYEUAAAAIVQQuoA3MhIlSbJyvUVYiu+5DR+sBAABAaCJwAW1g4hN8oauWZfEMAAAANIPABbSRmdhg8Yzd22V3b3esFgAAAIQmAhfQRqbvAGnQUH/bvrfUwWoAAAAQighcwGkwl07zb9u1H8oWH3GwGgAAAIQaAhdwGszYC6WU7r6GxyO76h/OFgQAAICQQuACToNxu2UmXuFv21V/l/VUO1cQAAAAQgqBCzhN5uKpkjva1zh2RHbdx84WBAAAgJBB4AJOk0nuKjP+Yn/bvvuWrLUOVgQAAIBQQeAC2oG57Mr6xu7tUs4254oBAABAyCBwAe3A9BsoZZ3lb9t333KwGgAAAIQKAhfQTlyT6me57PqPZQ8fcrAaAAAAhAICF9Bezh4v9ezt2/Z6ZVe86Ww9AAAAcByBC2gnxhUlM/kqf9t++E/Z0hIHKwIAAIDTCFxAOzIXfl1KTPY1KitkP+BFyAAAAJGMwAW0IxMbJzPxcn/bvvsWL0IGAACIYAQuoJ2Zy6ZJbrevcaxI9rMPnC0IAAAAjiFwAe3MdOkmc/5l/rb95+u8CBkAACBCEbiAIGi4eIby9kib1jlXDAAAABxD4AKCwKSfIY0+z9/2/v2vzHIBAABEIAIXECSuy/+tvrHzKyn7C+eKAQAAgCMIXECQmEFDpSEj/W3v2391sBoAAAA4gcAFBJHrimvrG19+Lrt7u3PFAAAAoMMRuIBgGjZayjzT3/T+/VUHiwEAAEBHI3ABQWSMaTzLteFT2bw9zhUEAACADkXgAoJt9HlS+hn+puVZLgAAgIhB4AKCzLhcMg1mueyaD2UP7HOwIgAAAHQUAhfQAcx5X5N69/E1rFd22SvOFgQAAIAOQeACOoBxRclMn+lv288+YJYLAAAgAhC4gA7CLBcAAEDkIXABHYRZLgAAgMhD4AI6kDnva1Jag1mupX9xtiAAAAAEFYEL6EC+Wa5v+9t2zQey+3kvFwAAQLgicAEdzJx7Uf17uayV940XnS0IAAAAQUPgAjqYcUXJdfX36jv+tVp2V7ZzBQEAACBoCFyAE845X+o/2N/0vv5HB4sBAABAsBC4AAcYY+SaMbu+Y+tG2a0bnSsIAAAAQUHgApwy7GxpyEh/0/v6H2Wtda4eAAAAtDu30wWEkpKSEt19990qLi5WRkaGnnjiiYDHrlq1SsuXL9e+ffvkdruVlZWlGTNmaMiQIR1XMDo1Y4xc35ot7/x7fR27sqXPP5XGXOBsYQAAAGg3zHA18MILL+j48eOndNzChQu1d+9ejRw5UoMHD9amTZv04IMPas2aNR1QKcKFGTRUGn2ev+197Q+yHo+DFQEAAKA9Ebhqbd68WatWrdKkSZNaPG7Lli1atmyZkpOT9ctf/lL33nuv5s6dq4ceekgul0uLFi1SSUlJB1WNcOCacb1kav8qHsqT/WC5swUBAACg3RC4JFVVVemZZ55R3759deWVV7Z47FtvvSVJmjFjhtLT0/39WVlZmjx5ssrKyrRy5cqg1ovwYjL6yXxtsr9t3/qzbFmpgxUBAACgvRC4JP31r3/VwYMHdfPNNysqKirgcVVVVdqyZYskacKECU321/WtX78+OIUibJlvfleKjfM1Sopll7/qbEEAAABoFxEfuHJzc7V06VJNnDhRw4cPb/HYvLw8VVdXq0uXLurRo0eT/QMGDPBfE2gN07WbzNQZ/rZd8ZZsUYGDFQEAAKA9RHTg8nq9evrpp5WQkKBZs2ad9PjCwkJJajZsSVJcXJwSExNVWlqq8vLydq0V4c9MuVrq2t3XqK6S5WXIAAAAnV5ELwu/fPly7dixQ7fffruSk5NPenxFRYUkKSYmJuAxsbGxKi0tVUVFheLj41u83pw5c5r0xcTEaP78+ZKknj17nrQmdAy32/dXJTU1Najfp2zW/1PxwkclSXb1++r6re8pZsiIoH5PnJ6OGhvonBgfCISxgUAYG+EnYme4CgsL9ec//1nDhw/XxIkTT+mcupfSGmNOegzQFvGXTZM780x/+/jiBbJer4MVAQAA4HRE7AzX4sWL5fF4dPPNN5/yOXUzVpWVlQGPqaqqkuS7vfBkFixY0OL+wsJCAlyIqPstU0FB8J+rstd+X/rlTyRJ1du36tBbr8h1QcuvK4BzOnJsoPNhfCAQxgYCYWyEHmNMo9XJWytiA9eGDRuUmJioxYsXN+qvrq6W5As78+bNkyTdd999iouL89/id/jw4WavWVFRodLSUiUmJp70dkIgEJM1QmbcRbLrPpIk2SV/kB1zvkxcgsOVAQAAoLUiNnBJUmlpqb788stm91VVVfn31dTUSJIyMjIUHR2t4uJiHT58uMniGbt27ZIk9evXL4hVIxKYf7tRduMaqbpKOnZEdtlfZa65wemyAAAA0EoRG7heeeWVZvsPHTqkH/7wh8rIyNATTzzRaF9MTIxGjBihzz//XKtXr9a0adMa7V+9erUkaezYsUGpGZHD9OglM3WG7NI/S5LsO3+TvXCSTFpfhysDAABAa0TsohltVReylixZovz8fH9/dna2VqxYofj4eF122WVOlYcwYr5xjdS9dqXKGo+8Lz/NM30AAACdTMTOcLXVqFGjdMUVV+jtt9/Wvffeq5EjR6qmpkabNm2S1+vVnXfeqaSkJKfLRBgwsbFyXXezvL/1vSZAWzfKrvtY5tyLnC0MAAAAp4zA1QY33nijMjMztXz5cm3evFlRUVEaMWKErrnmGg0dOtTp8hBOxpwvjRgjbdkgSbKvLJYdOYYFNAAAADoJY7lHKWTl5+dzC1mIcHKJVnsoT94H75Q8vhU0zeSr5Lrupg6vA81j+V60hPGBQBgbCISxEXpOd1l4nuECQpzplSFz+TX+tn33Ldl9uxysCAAAAKeKwAV0AuYb10ipab6G1yvvHxbKemucLQoAAAAnReACOgETEyvXd2+t79iVLfveMucKAgAAwCkhcAGdhBkxRmb8Jf62feNF2cOHHKwIAAAAJ0PgAjoRM/NmKSnZ16iskPel37KwCgAAQAgjcAGdiEnuKnPdzfUdm9fJrvnAuYIAAADQIgIX0MmYCROl4Wf72/bPz8gWH3WqHAAAALSAwAV0MsYYuWbdLsXE+jpKirm1EAAAIEQRuIBOyKSmyVxzQ33Hhk9k137oXEEAAABoFoEL6KTMxCukISP9bfvy07LHjjhXEAAAAJogcAGdlHG55LrhTik2ztdRelzeFxdxayEAAEAIIXABnZhJTZP5txvrO/71meyn7zlWDwAAABojcAGdnLn4G9Kw0f62ffl3sgUHHKwIAAAAdQhcQCdnXC65bvwPKSHR11FZLu+zj8vW1DhbGAAAAAhcQDgw3VNlZt1e37Fjq+zy15wrCAAAAJIIXEDYcJ37Nd9LkWvZt/4ku2u7cwUBAACAwAWEE/Od/yd1T/U1amrkfeaXsuVlzhYFAAAQwQhcQBgxCYly3TRHMrV/tQsOyP5xIUvFAwAAOITABYQZk3WWzPSZ/rZd+6HsR+84WBEAAEDkInABYchMv04aMtLftn/6nez+XAcrAgAAiEwELiAMGVeUXDfPkZK6+Dqqq+R9+heyFeXOFgYAABBhCFxAmDIpPeT697vrO/L3yv7h/3ieCwAAoAMRuIAwZkaOk/nGNf62Xfuh7HvLHKwIAAAgshC4gDBnrp7V+Hmuv/5edudXDlYEAAAQOQhcQJgzUVFy3XKPlNLd11FTI+9vfy5bfMTZwgAAACIAgQuIAKZLN7n+371SVJSv4+hheZ/6uayn2tnCAAAAwhyBC4gQZvBwmWv/vb5jx5eyf37GuYIAAAAiAIELiCDmsukyF0zyt+2q5fK+/3cHKwIAAAhvBC4gghhjZGbdJg3I8vfZP/9ONnuLg1UBAACELwIXEGFMdIxct/9Y6tpgEY2nHpU9lOdsYQAAAGGIwAVEIJPSwxe63NG+jpLj8v76YdnSEmcLAwAACDMELiBCmYFDZG78j/qOg/t9M12sXAgAANBuCFxABHONv0Tmqu/Wd2zbLPvHRbLWOlcUAABAGCFwARHOTJspM2Giv20/eVf2zT85VxAAAEAYIXABEc4YI3P9ndKZw/19dumf5V213MGqAAAAwgOBC4BMdLRcd8yV0s/w99mXfiv7+WoHqwIAAOj8CFwAJEkmMVmuu+ZJKT18HdYr7zOPyW7/0tG6AAAAOjMCFwA/0yNVrrselOITfR3VVfL+5n9kc3c6WxgAAEAnReAC0Ijpm+m7vbDuHV3lZfI+8aBs/j5nCwMAAOiECFwAmjBDRsh1631SVJSvo6RY3gU/lS086GxhAAAAnQyBC0CzzOhzZb5/t2SMr+PoYV/oKip0tC4AAIDOhMAFICDX+EtkvndbfUfBAXl/db/s0cPOFQUAANCJELgAtMh1yTdkrv1+fcehPF/oOnbEuaIAAAA6CQIXgJNyTfmWzIwb6jsO7K+d6SpyrigAAIBOgMAF4JS4Lr9G5upZ9R35e+X9xX0spAEAANACAheAU+aadp3MN79b31FwQN5f/Fj2AEvGAwAANMftdAFOqays1MaNG7V+/Xrt3LlTBQUF8nq9SktL0/jx4zV9+nTFxcU1e+6qVau0fPly7du3T263W1lZWZoxY4aGDBnSwZ8C6HiuK78tb0yM7KvP+zqOFMr7ix/L9Z//I3PGAEdrAwAACDURO8P10Ucf6bHHHtPKlStlrdXo0aM1dOhQHTp0SK+88op+/OMf69ixY03Oe+GFF7Rw4ULt3btXI0eO1ODBg7Vp0yY9+OCDWrNmjQOfBOh4rqkzZGbdXr9k/PFj8j72E9mdXzlbGAAAQIiJ2Bkut9utKVOmaNq0aUpPT/f3HzlyRPPnz9euXbv0/PPP66677vLv27Jli5YtW6bk5GQ98sgj/vOys7M1b948LVq0SMOHD1dSUlKHfx6go7ku+Ya8sXGyzz0heb1SWam8jz8g1x1zZYaNdro8AACAkBCxM1yXXHKJbr755kZhS5K6deumm266SZK0Zs0aeTwe/7633npLkjRjxoxG52VlZWny5MkqKyvTypUrO6B6IDS4JkyU69b7JHft724qK+T99f/IbmS2FwAAQIrgwNWS/v37S5Kqq6t1/PhxSVJVVZW2bNkiSZowYUKTc+r61q9f30FVAqHBnDNBrjt/KsXE+jo81fIu/Jm877/tbGEAAAAhgMDVjIMHfctcR0VF+W8PzMvLU3V1tbp06aIePXo0OWfAAN9iAbm5uR1XKBAizPBz5PrPh6T4BF+H9cq+9Ft5X/m9rLfG2eIAAAAcROBqxttv+34zf/bZZys6OlqSVFhYKEnNhi1JiouLU2JiokpLS1VeXt4xhQIhxAweLtd/PSql1P8dse/8Td6nfi5bWeFgZQAAAM6J2EUzAtmwYYNWrlypqKgozZw5099fUeH7B2NMTEzAc2NjY1VaWqqKigrFx8ef9HvNmTOnSV9MTIzmz58vSerZs2dry0eQuGufUUpNTXW4khCXmqqaXz2rI//7X/LkZPv6/rVariceVLe5v1RUt+Z/YdGZMTbQEsYHAmFsIBDGRvhhhquBffv26Te/+Y2stZo9e7YyMzP9+6y1kiRTtwx2M+qOASJZVI9e6v6/Tyl23AX+Ps+OrTp8782qzt3pYGUAAAAdjxmuWocPH9bPfvYzlZaWavr06briiisa7a+bsaqsrAx4jaqqKkkK+MLkEy1YsKDF/YWFhYS4EFH3W6aCggKHK+k87A/+SyZ5sezKZZIkb8EBHb7vFrlu+S+ZEWMdrq79MDbQEsYHAmFsIBDGRugxxjRZ2bw1mOGSVFxcrEceeUSFhYWaOHGiZs+e3eSYutv7Dh8+3Ow1KioqVFpaqsTExFO6nRAId8YVJdd3/5/Mt39Q/4Lk8jJ5f/0/8i79i6zX62yBAAAAHSDiA1d5ebkeffRR7d+/X+edd55uvfXWZm8bzMjIUHR0tIqLi5sNXbt27ZIk9evXL+g1A52Ja9KVct3+k/pl462V/dtL8i76mWxZibPFAQAABFlEB67q6mr94he/0M6dOzV69Gjdfffdcrma/58kJiZGI0aMkCStXr26yf66vrFjw+dWKaC9mLPHy/XjX0ipafWdG9fI+78/kt2327G6AAAAgi1iA5fX69WTTz6pL774QsOGDdM999zjXxUmkGnTpkmSlixZovz8fH9/dna2VqxYofj4eF122WVBrRvorEzfAXLdv0AadW5956F8eR+9R97PVjlXGAAAQBBF7KIZy5cv15o1ayRJycnJWrx4cbPHzZ49W126dJEkjRo1SldccYXefvtt3XvvvRo5cqRqamq0adMmeb1e3Xnnnf4XJQNoyiQkyXXHXNm3/yr75suStVJVleziX8mbs03m2u/LuKOdLhMAAKDdRGzgKimpf3akLng159prr/UHLkm68cYblZmZqeXLl2vz5s2KiorSiBEjdM0112jo0KFBrRkIB8blkpk+UzZzsLzP/EqqfY7LvrdUdsdWuX7wI5m0vg5XCQAA0D6MZd3xkJWfn8+y8CGCJVqDwxYckPe386U9OfWdMbEyM2+W+dqUFt97FyoYG2gJ4wOBMDYQCGMj9LAsPIBOy6SmyfXfP5e5+Bv1nVWVsn9cKO9Tj8qWFDtXHAAAQDsgcAFwlImJlWv27XLd8RMpKbl+x+er5X3oP2S//JdjtQEAAJwuAheAkGDOniDXg7+Whp9d33m0SN7HH5D3r8/KVlU6VhsAAEBbEbgAhAyT0kOuu+bJXPvvUoPXNNh/viHvQ3fJZn/hYHUAAACtR+ACEFKMyyXXlKvl+vFjUvoZ9TsO5cn7yx/L+9JvZSvKnCsQAACgFQhcAEKS6TdQrrkLZKZcLZn6H1X2/bflffCHspvXO1ccAADAKSJwAQhZJjZWrmv/Xa4f/0Lq079+R1GhvL9+SN7fP85KhgAAIKQRuACEPDMgS677F8hc+R0pqsGzXatXyvvT2+T94B+y3hoHKwQAAGgegQtAp2Dc0XJ98zty/fRxKfPM+h0lx33v7frZf8nu/Mq5AgEAAJpB4ALQqZg+/eX68S98KxnGxtXvyN0h7/x75X3+SdniI84VCAAA0ACBC0CnY1xRvpUMH35K5rxLGu2zH78r7/23ybviTVmPx6EKAQAAfAhcADot062HXD/4kVz/9TOpb2b9jvIy2b8slnfenbIbPpG11rEaAQBAZCNwAej0TNYIue5/XOY7t0gJifU7Du6X96n58j76X7LZW5wrEAAARCwCF4CwYKKi5LpsulyP/Fbma1MavbtLu7Ll/eVPVPObh2X35zpXJAAAiDgELgBhxSR3lev6H8r14K+l0ec13rlprbwP/Ye8zz0pW3DAmQIBAEBEcZ/8EADofEyffor64f2y27+U97Xnpbol462V/eRd2dUrZc6/VOaKa2V6ZThaKwAACF/McAEIa+bM4XL998/luv0nUlrf+h1eb+2KhrfL+/vHZQ/sd65IAAAQtpjhAhD2jDHSORPkGnWub3br7b9KhQd9O61XdvVK2c9WyZz7NZlp18pk9HO2YAAAEDYIXAAihomKkvnaFNnzL5P9bJXssr9Idc9yWa/smlWya1ZJo86Va8q3pKyzfGENAACgjQhcACKOcbtlLpwkO2Gi7JoPZJe9Ih1scEvhprXyblor9R8sM+VqmbEXykRFOVcwAADotHiGC0DEMlFRcp1/qVz/838yN/9IOvFWwtwdss88Ju9PbpH3nb/Jlpc5UygAAOi0mOECEPGMK0pm/CWy510sfbFB3n++IW3dWH9AUYHsK7+X/dvLMhMukZl4uUzfAU6VCwAAOhECFwDUMsZII8YqasRY2T05su/8TXbtB1JNje+AynLZVctlVy2XBg+TmXiF7NRvykTHOFo3AAAIXcZaa50uAs3Lz88X/3lCQ2pqqiSpoKDA4UrQ0WxRoex7S2U//IdUVtpkv6trN8V/fbrKx1wk0yvdgQoRyvjZgUAYGwiEsRF6jDFKT2/7/8cTuEIYgSt08MMPtrJSdt2HsivflnJ3NH9Q1lkyF0zyLbIRF9+xBSIk8bMDgTA2EAhjI/ScbuDilkIAOAUmNlbmwq9LF35ddtd22fffll37oVRdVX9Q9hey2V/I/ul3vtB14STpTJaWBwAgkhG4AKCVzIAzZQbcJXvdvytx42cqe+dN1ezbXX9AZYXvBcufvCulpsmcd7HvixcqAwAQcbilMIRxS2HoYHofgaSmpspaq4I1H/tC1poPpfKmz3pJkvr09wWvc78mk5rWsYXCEfzsQCCMDQTC2Ag93FIIAA4zxsgMHCIzcIjsdTfJfr5a9pP3pK3/khr+0mR/ruzrf5R9/Y/SgCyZ874mc875Mj16OVY7AAAILgIXALQjExMrM/4SafwlvhUO133om/U6caGNXdmyu7Jl//J7qd8gmXMmyIw5X0o/g2e+AAAIIwQuAAgS072nzJRvSVO+JXtgf334yt/b+MA9O2X37JT920tS7z6+8HX2eGnAmTKuKGeKBwAA7YLABQAdwKT1kZn+bdlpM6X9u2XXfCj7+afSgf2NDzy4X3b5a7LLX5OSusiMGCONGCszYoxMYrIzxQMAgDYjcAFABzLGSH0HyPQdIM24XjZ/r+yGT2U/X930tsOSYtnV70ur35c1LmnQUJlR42SGnyOdMUDG5XLkMwAAgFNH4AIAB5n0M2SmnSFNu0728CHfghv/+kza8aVUU1N/oPVKO76U3fGl7JI/SInJ0tCRMkNHywwfLaWm8+wXAAAhiMAFACHC9Ogl8/VvSl//pmxZqbT1X7Kb1sluXicdP9b44NLj0vpPZNd/IitJ3VNlho2Who2WGTpKpms3Jz4CAAA4AYELAEKQSUiUxl4oM/ZCWa/Xt7DGpnWyX2yQdm+XvN7GJxQVyH68Qvp4hS+A9ekvk3WWNHi4zODhMt17OvExAACIeAQuAAhxxuWSMs+UyTxT+uZ3fLNf27+Q3bpRdutGKW9P05P258ruz5VWvu0LYD16yQwaJp05TGbwcCmjH8+AAQDQAQhcANDJmIREafR5MqPPkyTZo0WyX22Stm6U/WqjVFTY9KTDh2QPH5LWrPIFsPhE3yIcg4fJDBoq9R8sE5/QoZ8DAIBIQOACgE7OpHSXmTBRmjBR1lrpUL4vgO3YKrvjS6nwYNOTykulLetlt6z3BTBjfO8AyzyzdjZtsNRvoEx0TAd/GgAAwguBCwDCiDFG6p0h0ztDuuQbkiR79HBt+Noqu/1Lae8u36qHDVkrHdgne2CftHqlL4RFRfmeBasLYf0HSen9ZKKjO/xzAQDQWRG4ACDMmZQe0riLZMZdJEmyFWVSTrZvifkdW33v/yorbXpiTY20J0d2T470wT/qQ1haX5kzBvjeJ3bGAN87wZK7duhnAgCgsyBwAUCEMXEJ0vCzZYafLUm+VRALDsjuypZ2b5fdvV3amyNVVTU9uaamfkEOve8LYZKU0r02gGVKfTJl+vSTevdlNgwAEPEIXAAQ4YzLVX8b4oSJkiRbUyPl7fGFr7oQtn+PVONp/iJHi6SjRbJb1vvO911Y6pUuZZwhk97P92dGPymtD8+GAQAiBoELANCEiYry3Sp4xgDpa1MkSdbj8T3ntXeXtG+X78+9u6SS4uYvYr3Swf3Swf2yn6/2dUm+IJaa5gtgvTOkXhkyvftIvTOkrt18z6EBABAmCFwAgFNi3G6pb6ZM30xJl0qSb1XEY0XS3t2y+3wBzObtkQ7sDzwbZr3SoTzpUJ7/lkT/rYmx8VLvdJleGb4A1qt25i2tj0xicnA/IAAAQUDgAgC0mTFGSukhpfSQGTnW3289HqnggO+2xPw9Ut5eXxA7uF/yBAhiklRZXr9QR9216jYSk323PvZMk3r2lnr2kunZ27fdPdU3KwcAQIghcAEA2p1xu6X0vlJ6Xxld4O+3NTVSQb4viB3YLx3Mkz2UJx3Mk44fa/mipcelnG2yOdvqr1e34XJJ3XpKPXvXh7CG29yqCABwCIGrDaqqqvTGG2/o448/VmFhoZKSkjR69GjNnDlTPXr0cLo8AAhZpnZZeaX11Ynxx5aVSAfzZQ/u991yeDBP9qDv1kOVl7V8Ya9XOnxIOnxIdtvm+mvWbbijpe49pW49Zbr3lLqlSt17ynRPre9PSGzPjwoAgCQCV6tVVVXp4Ycf1rZt29StWzeNGzdOBQUFev/997VhwwY98sgjSktLc7pMAOh0TEKSNOBMmQFnNuq31krHj0oHamfDCg9KhQdlCw9KhYd8z5CdjKdaOpQvHcqvD2FSo23FxftmyRoGsZQeMindfcved+0uJSb7VnUEAOAUEbha6fXXX9e2bduUlZWl+++/X3FxcZKkpUuX6g9/+IOeeuopPfTQQw5XCQDhwxgjdekmdekmk3VWk/22qlI6XFAfwg43CGOFB323Ip6KinIpf6+UvzdwKItyS127+b5SuvvCWFdfIDNdGwSzpGRuYQQASCJwtYrH49Hy5cslSTfddJM/bEnS9OnTtWrVKm3dulU5OTkaOHCgU2UCQEQxMbENnhdrypaX+W43PFIoe7hAOlIoFRXKHimUimrbLS3k0VCNx3dOUYHv2g2/T8Pj3G6pa3cd7tlLrpTu8sbGS11SpOSuUnKKTJcUqUtXX19CEuEMAMIYgasVvvrqK5WWlqp3794aMGBAk/3jx49Xbm6u1q1bR+ACgBBh4hOkvpm+Je2b2W+9XqnkmFTUMIgV+gLa0cO+lzofK5Kqqk79m3o80uFDqj58qNndjWfNomqDmC+AmeSU+jCW3NUXzpJTattdZNzRp14HAMBxBK5WyM3NlaRmw5Ykf8iqOw4AEPqMy+W/ZVGZZzYfyqz1LdxxrEg6WiRbF8KO1raPFUnHjvja1a0IZpJUU+O/jnRCGGumrfhEKSlZSuoiJXWRSUqWErv4+0xtv/+YxGTfqpEAAEfwE7gVCgsLJSngSoTdu3dvdBwAIDwYY6SERN9X+hnNhjKpNpiVlfpnxZK91fIeLVJJfp5UfFT2+FGp+KhvCfzjx3xhq7XKS31fBQd83/PEGpo7Jz7B9x6zhiGtNow1aiclSwnJUmKSFBPLrY4A0A4IXK1QUVEhSYqNjW12f90zXXXHncycOXOa9MXExGj+/PmSpJ49e7alTASBu/a3w6mpqQ5XglDD2EDzfHdC1I2PxGaeEbNer2zpcXmPHZH3aJG8x46opvZP74l/HjviexatrcrLfF+FB33f+8RamjvHHS1XUrJMUrJciclyJXeRSexS35dUu53o227UF9P8/0+iHj87EAhjI/wQuFrB2mb/L+mU9wMAUMe4XDLJXeVK7up7xuwkbGWFao4ekT1+VN7jx+QtPiZvcf22PV7XLpb3+FF5i4+1/vbGhjy+2TkdLVKr5+FiYuRK7FIb0pKbCWl1fybJlZgsk5AoV2KSTGKyTFw8S+8DCCsErlaIj4+XJFVWVja7v66/4eqFLVmwYEGL+wsLCwlxIaLut0wFBQUOV4JQw9hAS9p9fLjcUteevq8WGEkua6WqSqnkuFRSLJUWyx4v9rVLi319JcdlS4p9S+eX1Pa1ZnGQQKqq5K0qlPdIG26xN0aKS/Ddvhmf4PtKSPItfhKfIMUnSQn126bBMb4/E6XomJC/HZKfHQiEsRF6jDFKT09v8/kErlaou8Xv8OHDze4vKipqdBwAAE4xxkixcb6vHr5/wJ1KBLHV1VJZie+r1PdlG7bLAvfJU336hVtb/5xaw+5AhzfXGeWuD1/xif5tXzir7UvwbZu6/XEJUnx87Z8JPMMGoN0QuFqhf//+kqRdu3Y1uz8nJ6fRcQAAdDYmOrr+5c51fad4rq2qrA1fpb5Zs7ITglltOLNltfvLy3zBqqzUNxvXXmo89TN2DetrruZA1zCu+gAWF18byuJl6gLZiQEtLkGmUbv+XBMV1X6fDUCnQ+BqhaFDhyohIUEHDx7Url27miwP/9lnn0mSxowZ40R5AAA4ysTESjGxUkr9ar6nHNY8HqmizBe+GgQxW14mlZdIZbULf5SX1PbVHVtaf07NKb7A+pQK8vquW3ZqM20t7ouJrQ9pteHtSNcUmYREeeU6IcC1EOrc0cy6AZ0QgasV3G63vvGNb2jJkiV69tlnNXfuXP/zWkuXLlVubq6GDh2qwYMHO1wpAACdi3G7/cvWN+o/xfOttb5FQgIFtvKy2tDm27blZb6Zt7rwVlHuC3xtWar/ZKoqfV/Hjvi7WprPCxjcoty1M2cNvmLrQtqJ/XGN+2MbzNbVfjHzBnQMAlcrzZgxQ5s3b9a2bdt01113aejQoSosLNT27duVnJys22+/3ekSAQCIOMYY30xSTKyU0r2+vxXXsNb6nkMrL/OFr/LaEFZeJntCu26/rWjY12C7PW+RrFPj8d2KWXq8ae2BPlNL14uOaRTO/EEsNr5psKsNbU2DXd12nIyLAAc0h8DVSjExMXrwwQf1+uuv66OPPtLatWuVmJioSy65RDNnzmTBDAAAOiljjC+ERMdIXVIa72vltWxNTeMAdkJAS4pyyVtWqrLDhY1DXUX5CYGv3Hd7YzBUV/m+jh9rXHtLn6ul68XENhPgEk6YZYs/IdzFNR/uYuJ4PQDCBoGrDWJiYjRz5kzNnDnT6VIAAEAIMlFRUmKS76thf+2fibVLf1ecZOlvW7e8/4mzbpV14a35L1tZUd+ubLDvdN7NdjJ1t06e+Bla+nwtXa9hQGsQzEzDdqxvds1/C2VsM/tqjzfu6NP9hECbELgAAABCVKPl/dW98b42XM96PFKAMGZPDGct9de1Pe24UMmJKmu/T+MJuLYHuCh3gzAW1yjMmdi4+oAXG9cgqDWYhTsh3CmW5+BwaghcAAAAEcK43ZK76cyb1NYAV910hq020NXPwJU1DXYNjmsU4oKxaEmdGk/9O+ZO/BwtfcaWrumObjbA+UJc8wGu0b4TApxiY0/3UyIEEbgAAADQJsYdLSVFN1ldUmrDc291i5Y0N5N2YlDz91fU3z5Z2SDE1R1TFcRbKCVfvSXVTd75JrU9xB2IiZUrPkHemNjGM3KxDW+nPGFGrqV9MbE8D+cwAhcAAAAc12jRkuSuTfe34Zq2pqY2fDUIaZUVvqBWUVa/r1Fgq6i/jfLEABfs2yglqapS3gCrXLYpxBkjxcSdMJsW1/h2yUbPxZ083CkmhnfCtQKBCwAAAGHJREVJCYm+rxP3tfGa1lMtVVaeEODKaxcyaS7A1c3EBQhwFeWSN0grUUqStfXPw524q6XTWrqmcdUHuNgTbqU84dbKRrdS+m+jbLg/wdcXxi/2JnABAAAAp8i4o33PbrXXc3DW+mbNakNat4R42YpyHT2Q1/Itkw1Xo2xmBi9orxOQfNeue2n4ibtaOq2la7pc9UHshBk3Exsnc/6lMiPHnW7ljiBwAQAAAA7x3UoZ7ftK7qLo2lcGmO69fX+24Zq+1wlUnRDSfDNtqgxwy2TDVwo02VfR7AxZu/J6pfJS39eJn0eSBg1t86yk0whcAAAAQBjxvU4gNuCqh20KcV6v7z1rpxTgavcFDHC1fwZ4Vq1ZsXFtqDo0ELgAAAAAtMi4XPULbHTt1nhfG69pvTW+0FVxwq2RlRW1z8PVBzSTOfj0P4RDCFwAAAAAOpxxRdUumpHQdJ8D9QQLi/IDAAAAQJAQuAAAAAAgSAhcAAAAABAkBC4AAAAACBICFwAAAAAECYELAAAAAIKEwAUAAAAAQULgAgAAAIAgIXABAAAAQJAQuAAAAAAgSAhcAAAAABAkBC4AAAAACBICFwAAAAAECYELAAAAAIKEwAUAAAAAQULgAgAAAIAgIXABAAAAQJAQuAAAAAAgSAhcAAAAABAkBC4AAAAACBK30wUgMGOM0yXgBPw3QSCMDbSE8YFAGBsIhLEROk73v4Wx1tp2qgUAAAAA0AC3FAIAAABAkBC4gFNw33336b777nO6DIQgxgZawvhAIIwNBMLYCD88wwWcgqqqKqdLQIhibKAljA8EwthAIIyN8MMMFwAAAAAECYELAAAAAIKEwAUAAAAAQULgAgAAAIAg4T1cAAAAABAkzHABAAAAQJAQuAAAAAAgSAhcAAAAABAkBC4AAAAACBICFwAAAAAECYELAAAAAIKEwAUAAAAAQeJ2ugAgXLz66qt65ZVXJEl33XWXLrzwQocrghP279+vtWvXatOmTcrPz9exY8eUmJioIUOGaNq0aRo2bJjTJSLIqqqq9MYbb+jjjz9WYWGhkpKSNHr0aM2cOVM9evRwujw4pLKyUhs3btT69eu1c+dOFRQUyOv1Ki0tTePHj9f06dMVFxfndJkIESUlJbr77rtVXFysjIwMPfHEE06XhNNA4ALaQV5enl5//XUZY8S7xCPbww8/rKKiIsXHx+vMM89UVlaW9u3bpzVr1mjt2rW6/vrrNW3aNKfLRJBUVVXp4Ycf1rZt29StWzeNGzdOBQUFev/997VhwwY98sgjSktLc7pMOOCjjz7S008/LUk644wzNHr0aJWXlys7O1uvvPKKPv74Y82bN09du3Z1uFKEghdeeEHHjx93ugy0EwIXcJqstXr66aeVkJCgM888U+vWrXO6JDiob9++mjVrliZMmCC3u/5H7DvvvKNnnnlGf/zjHzV69Gj17dvXwSoRLK+//rq2bdumrKws3X///f4Zi6VLl+oPf/iDnnrqKT300EMOVwknuN1uTZkyRdOmTVN6erq//8iRI5o/f7527dql559/XnfddZeDVSIUbN68WatWrdLXv/51rVixwuly0A54hgs4Te+++662bt2q66+/XomJiU6XA4fdf//9uuiiixqFLUmaPHmyRo8eLa/Xq08//dSh6hBMHo9Hy5cvlyTddNNNjW4Pmz59uvr376+tW7cqJyfHqRLhoEsuuUQ333xzo7AlSd26ddNNN90kSVqzZo08Ho8T5SFEVFVV6ZlnnlHfvn115ZVXOl0O2gmBCzgNR48e1UsvvaSRI0fqa1/7mtPlIMT1799fku832gg/X331lUpLS9W7d28NGDCgyf7x48dLErPgaKLuZ0N1dTW3kUW4v/71rzp48KBuvvlmRUVFOV0O2gmBCzgNzz77rKqqqnTzzTc7XQo6gYMHD0qSUlJSnC0EQZGbmytJzYYtSRo4cGCj44A6dT8boqKilJSU5HA1cEpubq6WLl2qiRMnavjw4U6Xg3ZE4ALaaP369Vq9erW+9a1vNblFBDjRgQMHtGHDBknSuHHjHK4GwVBYWChJAVci7N69e6PjgDpvv/22JOnss89WdHS0w9XACV6v1/88+KxZs5wuB+2MwAW0QUVFhRYvXqz09HRdddVVTpeDEFdTU6NFixapurpaF1xwgX+mA+GloqJCkhQbG9vs/rpnuuqOAyRpw4YNWrlypaKiojRz5kyny4FDli9frh07dmj27NlKTk52uhy0M1YpRET61a9+pb1797bqnB/+8IcaPHiwJOnll1/W4cOH9cADD/DbyDBzumOjOc8++6y++uor9e7dm9tPw9jJXgnBKyNwon379uk3v/mNrLWaPXu2MjMznS4JDigsLNSf//xnDR8+XBMnTnS6HAQBgQsRqaCgQHl5ea06p7KyUpK0Y8cO/eMf/9DFF1+sESNGBKM8OOh0xkZzXn31Vb3zzjvq2rWr5s6dy/MZYSw+Pl5S4PFQ18/LbSFJhw8f1s9+9jOVlpZq+vTpuuKKK5wuCQ5ZvHixPB4Pv5ALYwQuRKT58+e3+dwNGzbIWqs9e/Zo3rx5jfbt379fUv0/sidMmKBvfOMbp1MqOtjpjI0TLV++XK+88ooSEhI0d+5cXngb5nr27CnJ9w/p5hQVFTU6DpGruLhYjzzyiAoLCzVx4kTNnj3b6ZLgoA0bNigxMVGLFy9u1F9dXS3JNwNW9++N++67j1/adEIELqCNdu/eHXDf/v37tX//fm4PiWAffvihnnvuOcXGxuq+++5jLESAuqW9d+3a1ez+uvdv1R2HyFReXq5HH31U+/fv13nnnadbb71Vxhiny4LDSktL9eWXXza7r6qqyr+vpqamI8tCOyFwAa103XXX6brrrmt238KFC7Vq1SrddddduvDCCzu4MoSKDRs2aNGiRYqKitI999yjoUOHOl0SOsDQoUOVkJCggwcPateuXU2Wh//ss88kSWPGjHGiPISA6upq/eIXv9DOnTs1evRo3X333XK5WL8s0r3yyivN9h86dEg//OEPlZGRoSeeeKJji0K74m85ALSjr776SgsWLJAk3X333Ro9erTDFaGjuN1u/y3Ezz77bKPVCJcuXarc3FwNHTq0xQVWEL68Xq+efPJJffHFFxo2bJjuueceud383huIBPxNB4B29POf/1xVVVXq1auX1q5dq7Vr1zY5ZujQoZo0aZID1SHYZsyYoc2bN2vbtm266667NHToUBUWFmr79u1KTk7W7bff7nSJcMjy5cu1Zs0aSVJycnKT53XqzJ49W126dOnI0gAEGYELANpRaWmpJN+tIIcOHQp4HIErPMXExOjBBx/U66+/ro8++khr165VYmKiLrnkEs2cOZMFMyJYSUmJf7sueDXn2muvJXABYcZYXgwCAAAAAEHBM1wAAAAAECQELgAAAAAIEgIXAAAAAAQJgQsAAAAAgoTABQAAAABBQuACAAAAgCAhcAEAAABAkBC4AAAAACBICFwAAAAAECQELgAAAAAIEgIXAAAAAAQJgQsAAAAAgoTABQAAAABBQuACAAAAgCAhcAEAAABAkBC4AAAAACBICFwAALSDN954Q9ddd52+853vaMeOHc0es2HDBs2cOVPXXXedPvroow6uEADgBAIXAADt4KqrrtLIkSNVU1OjJ598UuXl5Y32HzlyRIsWLZK1VhdffLEuuugihyoFAHQkAhcAAO3AGKM777xTXbt21cGDB/XMM8/491lr9X//938qLi5WWlqabr75ZgcrBQB0JAIXAADtJCUlRbfffruMMfroo4/0/vvvS5L+9re/afPmzYqKitJdd92luLg4ZwsFAHQYAhcAAO3onHPO0bRp0yRJzz77rD744AP95S9/kSR95zvf0aBBg5wsDwDQwYy11jpdBAAA4cTj8ej+++9XTk6Ov2/06NH6yU9+ImOMg5UBADoaM1wAALQzt9ut22+/3d9OSEjQHXfcQdgCgAhE4AIAIAhWrFjh3y4vL9fu3budKwYA4BgCFwAA7Wz9+vVavny5JKl///6y1mrhwoU6evSos4UBADocgQsAgHZU974tSZo4caIeeughpaam6tixY1q4cKF4dBoAIguBCwCAduL1evV///d/On78uNLT0/Xv//7vSkhI0F133aWoqCht3LhRS5cudbpMAEAHInABANBO3nzzzWbft5WVlaV/+7d/kyT96U9/arR6IQAgvBG4AABoBzt27Gj0vq2BAwc22v+tb31LZ511ljwej5588klVVFQ4USYAoIMRuAAAOE3l5eV68sknVVNTo1GjRunKK69scozL5dIPf/hDJSUlKT8/X88++6wDlQIAOhovPgYAAACAIGGGCwAAAACChMAFAAAAAEFC4AIAAACAICFwAQAAAECQELgAAAAAIEgIXAAAAAAQJAQuAAAAAAgSAhcAAAAABAmBCwAAAACChMAFAAAAAEFC4AIAAACAICFwAQAAAECQELgAAAAAIEgIXAAAAAAQJAQuAAAAAAgSAhcAAAAABAmBCwAAAACChMAFAAAAAEFC4AIAAACAIPn/IYN0FbuhsykAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 960x720 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "x = np.linspace(-5, 5, 201)\n",
    "\n",
    "plt.plot(x, np.exp(-x) - x)\n",
    "plt.xlabel('x')\n",
    "plt.ylabel('y')\n",
    "\n",
    "plt.legend([\"$e^{-x} - x$\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "79126eb7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 0.5000, err_a = 1.000e+00\n"
     ]
    }
   ],
   "source": [
    "# First Step\n",
    "x0 = 0\n",
    "x1 = x0 - (np.exp(-x0) - x0)  / (-np.exp(-x0) - 1)\n",
    "print(\"x1 = {:.4f}, err_a = {:.3e}\".format(x1, abs((x1 - x0)/x1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "26f21435",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 0.5663, err_a = 1.171e-01\n"
     ]
    }
   ],
   "source": [
    "# Second step\n",
    "x0 = x1\n",
    "x1 = x0 - (np.exp(-x0) - x0)  / (-np.exp(-x0) - 1)\n",
    "print(\"x1 = {:.4f}, err_a = {:.3e}\".format(x1, abs((x1 - x0)/x1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "c1788e25",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 0.5671, err_a = 1.467e-03\n"
     ]
    }
   ],
   "source": [
    "# Thrid step\n",
    "x0 = x1\n",
    "x1 = x0 - (np.exp(-x0) - x0)  / (-np.exp(-x0) - 1)\n",
    "print(\"x1 = {:.4f}, err_a = {:.3e}\".format(x1, abs((x1 - x0)/x1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "537b6504",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 0.5671, err_a = 2.211e-07\n"
     ]
    }
   ],
   "source": [
    "# Fourth step\n",
    "x0 = x1\n",
    "x1 = x0 - (np.exp(-x0) - x0)  / (-np.exp(-x0) - 1)\n",
    "print(\"x1 = {:.4f}, err_a = {:.3e}\".format(x1, abs((x1 - x0)/x1)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "66589001",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      converged: True\n",
       "           flag: 'converged'\n",
       " function_calls: 41\n",
       "     iterations: 39\n",
       "           root: 0.5671432904109679"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Bisection Method\n",
    "from scipy.optimize import root_scalar\n",
    "\n",
    "f = lambda x: np.exp(-x) - x\n",
    "root_scalar(f, bracket=[0, 1], method='bisect')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cfbdc1b9",
   "metadata": {},
   "source": [
    "### Quadratic Convergence \n",
    "엄밀해를 $x_r$ 로 가정하자 ($f(x_r) = 0)$. $x_i$를 중심으로 한 Taylor Expansion식을 적용하면 다음과 같다.\n",
    "\n",
    "$$\n",
    "0 = f(x_{r}) = f(x_{i}) + f'(x_i)(x_{r} - x_i) + \\frac{f''(\\xi)}{2!} (x_{r} - x_i)^2.\n",
    "$$\n",
    "\n",
    "이 식에서 Newton-Raphson 식을 빼면 다음과 같다.\n",
    "\n",
    "$$\n",
    "\\begin{align}\n",
    "0 &= f(x_{i}) + f'(x_i)(x_{r} - x_i) + \\frac{f''(\\xi)}{2!} (x_{r} - x_i)^2 \\\\\n",
    "&- f(x_{i}) - f'(x_i)(x_{i+1} - x_i) \\\\\n",
    "&= f'(x_i)(x_r - x_{i+1}) + \\frac{f''(\\xi)}{2!} (x_{r} - x_i)^2\n",
    "\\end{align}\n",
    "$$\n",
    "\n",
    "엄밀해로 부터 오차를 다음과 같이 정의하면\n",
    "\n",
    "$$\n",
    "E_i = x_r - x_i,\n",
    "$$\n",
    "\n",
    "위 식을 아래와 같이 정리할 수 있다.\n",
    "\n",
    "$$\n",
    "0 = f'(x_i) E_{i+1} + \\frac{f''(\\xi)}{2!} E_i^2\n",
    "$$\n",
    "\n",
    "해가 수렴하는 경우 $x_i, \\eta \\rightarrow x_r$ 이므로,\n",
    "\n",
    "$$\n",
    "E_{i+1} \\approx \\frac{-f''(x_r)}{2f'(x_r)} E_i^2 = O(E_i^2).\n",
    "$$\n",
    "\n",
    "즉 현재의 오차는 이전 단계에서의 오차에 제곱에 비례한다. \n",
    "\n",
    "#### 예제\n",
    "위 예제에서 참값은 $x_r=0.5671432904109679$ 이다. 오차 $E_i$ 의 변화를 계산해보자."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "93ad8a06",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5671432904109679\n",
      "0.06714329041096789\n",
      "0.0008322872137497273\n",
      "1.253761057196101e-07\n",
      "1.1868284133242923e-12\n"
     ]
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : np.exp(-x) - x\n",
    "fp = lambda x: -np.exp(-x) - 1\n",
    "\n",
    "x0 = x = 0\n",
    "xr = 0.5671432904109679 \n",
    "\n",
    "Err = []\n",
    "itmax = 5\n",
    "for i in range(itmax):\n",
    "    # Newton-Raphson\n",
    "    xn = x - f(x) / fp(x)\n",
    "    \n",
    "    # Compute \n",
    "    Ei = xr - x\n",
    "    Err.append(Ei)\n",
    "    \n",
    "    x = xn\n",
    "    print(Ei)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "51f697cd",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.1809481283173079"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "-np.exp(-xr) / (2*fp(xr))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b37dfde7",
   "metadata": {},
   "source": [
    "위 오차식을 예제에 적용하면\n",
    "\n",
    "$$\n",
    "E_{i+1} \\approx \\frac{-f''(x_r)}{2f'(x_r)} E_i^2 = 0.18095 E_i^2.\n",
    "$$\n",
    "\n",
    "초기 오차 $E_0 = x_r - 0$ 이므로, 이 식으로 계산한 오차는 다음과 같다.\n",
    "\n",
    "- 매 계산마다 정확한 유효자릿수가 2개가 된다."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "5bfe9e19",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.06714329041096789 0.058202841070737574\n",
      "0.0008322872137497273 0.0008157626708729342\n",
      "1.253761057196101e-07 1.2534442801669386e-07\n",
      "1.1868284133242923e-12 2.8443834288658173e-15\n"
     ]
    }
   ],
   "source": [
    "# 실제 오차, 추정식\n",
    "for i in range(1, 5):\n",
    "    print(Err[i], 0.18095*Err[i-1]**2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6b74d052",
   "metadata": {},
   "source": [
    "### 예제 2 \n",
    "아래 함수의 근을 구하시오. 초기값은 0.5로 하자.\n",
    "\n",
    "$$\n",
    "f(x) = x^{10} - 1\n",
    "$$\n",
    "\n",
    "Newton-Raphson 법을 적용하면 다음과 같다.\n",
    "\n",
    "$$\n",
    "x_{i+1} = x_{i} - \\frac{f(x_i)}{f'(x_i)} = x_i - \\frac{x_i^{10} - 1}{10 x_i^9}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "79ea1fb2",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 51.6500, err = 5.115e+01\n",
      "x1 = 46.4850, err = -5.165e+00\n",
      "x1 = 41.8365, err = -4.648e+00\n",
      "x1 = 37.6529, err = -4.184e+00\n",
      "x1 = 33.8876, err = -3.765e+00\n",
      "x1 = 30.4988, err = -3.389e+00\n",
      "x1 = 27.4489, err = -3.050e+00\n",
      "x1 = 24.7040, err = -2.745e+00\n",
      "x1 = 22.2336, err = -2.470e+00\n",
      "x1 = 20.0103, err = -2.223e+00\n",
      "x1 = 18.0092, err = -2.001e+00\n",
      "x1 = 16.2083, err = -1.801e+00\n",
      "x1 = 14.5875, err = -1.621e+00\n",
      "x1 = 13.1287, err = -1.459e+00\n",
      "x1 = 11.8159, err = -1.313e+00\n",
      "x1 = 10.6343, err = -1.182e+00\n",
      "x1 = 9.5708, err = -1.063e+00\n",
      "x1 = 8.6138, err = -9.571e-01\n",
      "x1 = 7.7524, err = -8.614e-01\n",
      "x1 = 6.9771, err = -7.752e-01\n",
      "x1 = 6.2794, err = -6.977e-01\n",
      "x1 = 5.6515, err = -6.279e-01\n",
      "x1 = 5.0863, err = -5.651e-01\n",
      "x1 = 4.5777, err = -5.086e-01\n",
      "x1 = 4.1199, err = -4.578e-01\n",
      "x1 = 3.7079, err = -4.120e-01\n",
      "x1 = 3.3371, err = -3.708e-01\n",
      "x1 = 3.0034, err = -3.337e-01\n",
      "x1 = 2.7031, err = -3.003e-01\n",
      "x1 = 2.4328, err = -2.703e-01\n",
      "x1 = 2.1896, err = -2.432e-01\n",
      "x1 = 1.9707, err = -2.189e-01\n",
      "x1 = 1.7738, err = -1.968e-01\n",
      "x1 = 1.5970, err = -1.768e-01\n",
      "x1 = 1.4388, err = -1.582e-01\n",
      "x1 = 1.2987, err = -1.401e-01\n",
      "x1 = 1.1784, err = -1.204e-01\n",
      "x1 = 1.0833, err = -9.500e-02\n",
      "x1 = 1.0237, err = -5.969e-02\n",
      "x1 = 1.0023, err = -2.135e-02\n",
      "x1 = 1.0000, err = -2.292e-03\n",
      "x1 = 1.0000, err = -2.393e-05\n",
      "x1 = 1.0000, err = -2.578e-09\n",
      "x1 = 1.0000, err = 0.000e+00\n",
      "x1 = 1.0000, err = 0.000e+00\n"
     ]
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : x**10 - 1\n",
    "fp = lambda x: 10*x**9\n",
    "\n",
    "x0 = x = 0.5\n",
    "\n",
    "Err = []\n",
    "itmax = 45\n",
    "for i in range(itmax):\n",
    "    # Newton-Raphson\n",
    "    xn = x - f(x) / fp(x)\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "610497dc",
   "metadata": {},
   "source": [
    "### 예제 3\n",
    "아래 함수의 해를 적절한 초기값으로 부터 근을 구하시오.\n",
    "\n",
    "$$\n",
    "f(x) = x^2 + 10 \\sin(x)\n",
    "$$\n",
    "\n",
    "Newton-Raphson 법을 적용하면 다음과 같다.\n",
    "\n",
    "$$\n",
    "x_{i+1} = x_{i} - \\frac{f(x_i)}{f'(x_i)} = x_i - \\frac{x^2 + 10\\sin(x)}{2 x + 10\\cos(x)}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "348396d1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x7fe95bfb1210>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1wAAAKLCAYAAADiqYnnAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAABcSAAAXEgFnn9JSAACV10lEQVR4nOzdd3yV5f3/8fd1MkhIAshesg1D9gYHGxVwoWLrbmuto1bL11pr+6vaWke/fh11r7baahUtDoYoIOBir7AEGbJXWIGQfV+/Pw45J5EwAufkOuP1fDx8eH+uM/LOlQD55L7v6zLWWisAAAAAQMj5XAcAAAAAgFhFwwUAAAAAYULDBQAAAABhQsMFAAAAAGFCwwUAAAAAYULDBQAAAABhQsMFAAAAAGFCwwUAAAAAYULDBQAAAABhQsMFAAAAAGFCwwUAAAAAYULDBQAAAABhQsMFAAAAAGFCwwUAAAAAYULDBQAAAABhkug6AI5tx44dstY6+dh169aVJGVnZzv5+LGIOQ095jT0mNPQY05DjzkNPeY09JjT0HM1p8YYNWzY8JRfT8MVway1zhqushkQWsxp6DGnocechh5zGnrMaegxp6HHnIZetM0plxQCAAAAQJjQcAEAAABAmNBwAQAAAECY0HABAAAAQJjQcAEAAABAmNBwAQAAAECYsCw8AABAjDjV5bJLXxdty21HMuY09EI1p8aYUMQ5aTRcAAAAUay4uFh5eXkqLCw85R9E9+3bJ0nyPC+U0eIacxp6oZpTY4ySk5OVmpqqxMTwt0M0XAAAAFGquLhYBw4cUEpKimrVqiWf79TuFin9obO4uDiU8eIacxp6oZpTz/OUn5+vAwcOqGbNmmFvumi4AAAAolReXp5SUlKUlpZ2Wu9TeolVVV9qFcuY09AL1ZwmJCQE/szk5eUpIyPjtLMdD4tmAAAARKnCwkKlpKS4jgFEpWrVqqmwsDDsH4eGCwAAIApZa2WtPeXLCIF4l5CQEPhzFE78CQUAAACAMKHhAgAAAIAwoeECAAAAgDCh4QIAAACAMKHhAgAAAIAwoeECAAAAgDCh4QIAAEDcevfdd3XllVeqU6dOyszM1IUXXqjx48e7joUYkug6ACKbtZYd0gEAQMz66quvNHz4cP3hD39QzZo1NWXKFP3qV79SQkKCLr30UtfxEAM4w4Wj2D27dej9N7T3//1S3nMPu44DAAAQNs8++6xuueUWde3aVS1bttRtt92mQYMGaeLEic4yZWVl6bnnntPNN9+sHj16qEmTJmrVqtUJX5efn68nnnhC5557rlq1aqXu3btr7Nix2rZt22lnatKkifr06XPa71PWk08+qTPPPFNr1qw5rffZuXOnWrdurfvvvz9EyUKLM1w42t7dOvTvl/zHiUmyhQUyydXcZgIAAKgiBw4cUKNGjZx9/KefflqffvpppV6Tn5+vq6++WgsWLFCDBg00fPhwbdmyRe+++66mTZumjz/+WC1atAhP4FOwe/duvfTSSxo1apQyMzNP670aNGiga6+9Vm+88YZ++tOfqk2bNiFKGRqc4cLRWmbKpFT3HxcXSeu+dZsHAACgiowbN05ZWVm6/vrrnWXo0aOHfv3rX+uf//ynlixZclKvefbZZ7VgwQL16NFDX375pV566SVNnDhRf/zjH7Vnzx79z//8z2llmjVrlt59993Teo+ynn32WeXm5uqXv/xlSN7vtttuk+d5euKJJ0LyfqFEw4WjmMREJXfsGqjtqiXOsgAAAFSVTz/9VPfdd58ef/xxderU6ZTeY/PmzWrSpIl+9atfnXKOO+64Q/fcc4+GDRumevXqnfD5RUVF+sc//iFJeuSRR5SWlhZ47Be/+IXat2+vOXPmKCsr65QztWnTJmRnyPLy8vTee++pffv2Ovvss0Pyno0aNVL//v01ZcoU7d69OyTvGSo0XKhQcudegWO7cqnDJAAAACfvjjvuUJMmTfTMM88c9di8efPUqlUrde7cWd9//325xz766CPddtttevTRR3X11VdXUdrQmDdvng4cOKAWLVqoY8eORz0+cuRISdLUqVPLjX/33Xe688471b9/f7Vq1UqdOnXSsGHD9Mc//lE7d+4s99yK7uEqbS6vvPJK5eXl6ZFHHlHv3r3VsmVLnXPOOXr++edlrT0qz4QJE5STk6PRo0dX+Pkc72s4d+7cY34NL7/8chUVFWncuHEVvq8rNFyoUHKXYMOlTetkD+W4CwMAAHCS7rnnHiUmJuqVV17RwYMHA+Nr167VT37yE/l8Pr355pvlzta89dZb+vWvf62nnnoq6potSVq5cqUkVdhsSQqcrSt9niQtW7ZMF154oT744APVrl1bF154obp166aioiK9/vrrWrdu3Ul//MLCQl1zzTV66623dNZZZ6l///7asWOHHnnkEf31r3896vnTpk2TJPXr16/C9zve1/DGG2+s8GtY9v2mT59+0tmrAg0XKpTYrJV8tWr7C2ul1cvcBgIAADgJLVu21NVXX639+/frtddek+RfoOH666/XwYMH9dJLL6lr166B57/88sv63e9+p4ceekj9+vXTrl27tGvXLu3Zs8fRZ1B5W7dulaRjLvRROl76PEl6/fXXlZ+fr1deeUUTJ07UCy+8oDfffFMzZ87UzJkz1bp165P++AsXLpQxRl9++aXeeustvfXWWxo/frwSExP16quvKjc3t9zz58+fr6SkJHXo0KHC9zvW1/DHP/6xcnJyjvoalmrevLlq166tJUuWqKCg4KTzhxsNFypkjCl3lovLCgEAQLS4++67lZKSoldffVU7duzQjTfeqE2bNumxxx7T0KFDyz3373//u0pKSnTfffepW7dugf9GjBjhKH3lHT58WJKUmppa4ePVq/sXQyvb+JQ2lOecc85Rzz/rrLPUoEGDk/74Pp9P//u//6vatWsHxrp06aJBgwYpLy9PS5cGf47Mzs7Wrl271LRpU1WrduxVsI/1NfzrX/961NewrNatW6ugoKBSZ+jCjWXhcUzJnXspf5Z/SVIWzgAAIHpYa6W83BM/sfT5if4fCW1xcbginZzUNBljTvttGjdurOuvv16vvvqqhg0bpr1792rs2LG65pprjnru3LlzT/njPPfcc1q7dm25sdKmZu7cufrVr34lz/PKPf7LX/4y5MuWl94nday5q+g+qk6dOunzzz/XXXfdpbvuuktdunSRz3dq52LOPPPMCs+Ile4dtmvXrsBYdna2JKlWrVrHfc+Kvob33HOPrrvuOhUf5/u09H0j6QwlDReOqVrnnsFi9w7Z3Ttk6jV0FwgAAJycvFx5dx3dXBxLYRijVIbvmbel6ukhea9bbrlFr732mvbu3asxY8ac9rLoFZk5c6Zmz55d4WPff//9UYs6SNKYMWNC3nCVrkpYeqbrh/Ly8so9T/Ivoz5v3jxNnTpVU6dOVY0aNdStWzcNHTpUY8aMUXr6yX8djnUpY+nHK3t5X+k9WWWzHMsPv4a/+c1vTviajIyMch8nEtBw4ZgS6jWQGjaRdviv97XfZtFwAQCAiGet1UMPPRQ4s5OYGJ4fed9///2jxjZv3qy+ffvq6quv1t/+9rfjno0JlSZNmkiStm/fXuHjpeOlz5P8jcl7772n+fPna+rUqZo9e7a++uorzZo1S88995zGjx9/0svAV+asZGlDdOjQoeM+71S/hjk5OeU+TiTgHi4cl2nfJVis4j4uAAAQ+R588EFNnDhRw4YNU926dTVu3DitX7/edaywKV18Yvny5RU+vmyZf/Gz9u3blxs3xqh37976/e9/r4kTJ2rx4sW67LLLtHPnTj322GNhyVq3bl1J0v79+4/7vFP9Gh44cECSVKdOndPOGiqc4cJxmfZdZWdMliTZVUtlPU/mFK/vBQAAVSQ1zX953kkqPXtQFWdjjiv1xJeZncjLL7+s1157Td26ddOLL76oN998U3/605/0xBNP6IUXXghByMjTq1cv1ahRQ99//72WL19+1PLwkyZNkqTjLjYh+ZuUsWPH6sMPP9S3334blqx169ZV/fr1tWXLFuXl5VW40MexvoaPP/64Xn755eO+/9q1a5WSklKpVRbDjZ+ccXxtO0rmyLfJoRxpy/dO4wAAgBMzxshUT4++/05zwYyPP/5Yf/7zn9WiRQu98cYbSk1N1Q033KD69evr448/LrcPVSxJTk7WTTfdJEn6/e9/X+5erpdfflmrVq1S7969yy2l/uabb2rTpk1HvdeMGTMklb/8MNR69+6t4uLiCs/IHe9r+NFHH2nFihXHfN/vv/9e+/btU9euXY+7AmJVo+HCcZnq6VKL4I2dlssKAQBABJo9e7buvvtu1a5dW2+99VbgkrLU1FTdcccdstbq8ccfd5zy5EybNk2jRo0K/Cf5NxcuO1a6eXCpu+66S926ddOCBQt07rnn6tZbb9WoUaP0pz/9SWeccYaefPLJcs//17/+pX79+mnQoEH6+c9/rttuu03Dhw/XAw88oJSUFP36178O2+c3ZMgQSdI333xTbvxkvoaPPvroMd+3dAGTwYMHhyn5qaHhwgmZ9l0DxywPDwAAIs2aNWv0s5/9TD6fT2+88cZRiz1cd911atiwoaZNm6YFCxa4CVkJe/bs0eLFiwP/Sf5FJMqO/XDZ85SUFL333nu6++67lZqaqk8//VRbtmzRVVddpU8//VQtW7Ys9/zf/OY3+tGPfiRjjL7++mtNnTpVeXl5uvbaazVt2jT17NlT4XLxxRerRo0a+vDDDwNjJ/M1bNSokaZOnXrMr+EHH3ygpKQkjRkzJmzZT4WxFS3Mj4iwffv2CvdNqAr16tWT5N/V265eJu+J3/sfSE6W7+n/yCQlOckVzcrOKUKDOQ095jT0mNPQY079rLXas2eP6tSpc9qX4kXMPVwxhDk9vgceeECvvfaapkyZok6dOp3Ua443p9u2bVOfPn00cuRIvfTSSyf1fif7Z8gYc8yl708GZ7hwYq3aScnJ/uPCQml9eG6iBAAAQHy48847lZaWpueeey4k7/fSSy/J5/PpnnvuCcn7hVLMrFK4fv16ZWVlae3atfruu++0b98+JSUl6a233jru62bNmqUpU6Zoy5YtSkxMVGZmpkaPHq22bdse8zWrV6/W+PHjtWbNGhUXF6tp06a64IILNHDgwBB/VpHBJCVJZ50trThySnvlUpm2J/ebCAAAAOCH6tatq9tuu01PPvmk1qxZo8zMzFN+r507d+qtt97SNddcE/JNpUMhZhqu999/v9LX5L7xxhuaNGmSkpOT1blzZxUVFSkrK0tLly7V2LFj1bt376NeM2/ePD355JOy1qp9+/bKyMjQ8uXL9cILL2jjxo268cYbQ/UpRRTTvqtsacO1aol0+XVuAwEAACCq/frXvw7J4hwNGjTQunXrQpAoPGKm4crMzFSLFi3UunVrtW7dWrfccstxn798+XJNmjRJGRkZevjhhwPXZa5Zs0YPPvigXnjhBXXo0EHp6emB1xw6dEgvvPCCPM/T//zP/6hPnz6S/Bu3/fGPf9SkSZPUo0ePo/Y+iAWmfRcF7ib7fq3s4UP+FQwBAAAAHFPM3MN12WWXacyYMerRo4dq1ap1wudPmDBBkjR69OhyN8FlZmZq2LBhOnz4cGAfglKff/65Dh8+rJ49ewaaLUmqVauWrrvOf8Zn4sSJIfhsIlDTFlJ6Df+x9aTVFe9kDgAAACAoZhquyigsLAxstNa3b9+jHi8dW7hwYbnx0rqi13Tv3l1JSUlatmyZCgsLQx3ZOePzybTvEqhZHh4AAAA4sbhsuLZt26aioiLVqFEjsKFaWaX7FGzcuLHceOlu3K1atTrqNYmJiWrWrJmKioq0bdu2MKSOAOUaLjZABgAAAE4kLhuu7OxsSaqw2ZL8G8elpaUpNzdXeXl5kqTDhw8rNzdXklS7du0KX1c6Xvr+sabsGS7t2Cq7N773PwEAAABOJGYWzaiM/Px8SVJy6d5SFahWrZpyc3OVn5+v1NTUwGtKHzvWa8q+/4mMHTv2qLHk5GQ99thjkvzLZbpSurFc6eaSRwrtbthEJTu2SpLSt2xQ9bYdXMSLShXOKU4Lcxp6zGnoMaehx5z6WWu1b98+JSYmnvbGx6WvL51bnD7mNPRCPafWWvl8PtWrV++0/wwdT1ye4bLWv97e8Sa29DkoL7lrcKn8wqz5DpMAAAAAkS8uW+7U1FRJUkFBwTGfU7rwRUpKSrn/l76uevXqR72m9P3KPvd4nnzyyeM+np2d7azxK/2t4e7d5S8btC2CG0LnL56rXbt2hfU3ArHkWHOKU8echh5zGnrMaegxp0HWWhUWFiohIeG03qf0jEFxcXEoYkHMaTiEek6Li4tlrT3h7UDGmHKrmldWXJ7hKr1Ub8+ePRU+np+fr9zcXKWlpQWas+rVqwearL1791b4utJxl5cChl27TlJpg5WzX9q2yWkcAADiWXJy8knfygCgvIKCguPeYhQqcXmGq3HjxkpKSlJOTo727Nlz1OIZGzZskCQ1a9as3Hjz5s21atUqrV+/Xk2bNi33WHFxsTZt2qSkpCQ1btw4vJ+AQyYtQ2rWWtq4VpJ/eXjTpLnjVAAAxKfU1FQdOHBAkv9e8lM901V6RQ23VIQOcxp6oZrTkpISFRQUKD8/XzVr1gxFtOOKy4YrOTlZHTt21OLFizVnzhyNHDmy3ONz5syRJPXo0aPcePfu3bVq1SrNmTNH559/frnHFi1apKKiInXr1q1KOmWXTIcusqUN18ql0tBLHScCACA+JSYmqmbNmsrLy9OBAwdO+QdRn89/0ZPneaGMF9eY09AL1ZwaY5ScnKyaNWtWyaImcdlwSdLIkSO1ePFijR8/Xt27dw9cl7lmzRpNmzZNqampGjx4cLnXDBkyROPHj9eCBQs0d+5c9enTR5J04MAB/fvf/5YkjRo1qmo/EQdM+66yn/zXX6xZLltcLMMKPAAAOJGYmKiMjAxJp/6bf+6LCz3mNPRCNadVvf5AzPyUvGjRIv33v/8tN1ZcXKzf//73gfqKK65Q9+7dJUmdO3fWiBEjNHnyZN17773q1KmTSkpKlJWVJc/zdOeddyo9Pb3c+6Wnp+u2227TU089pSeffFIdOnRQRkaGli1bptzcXF100UXq1KlT+D9Z19q0l5KSpaJCqSBf2rBGOovl4QEAcO1Uf5AsfR0LYYUOcxp60TqnMdNw5eTk6Lvvvis3Zq0tN5aTk1Pu8ZtuukktWrTQlClTtGzZMiUkJKhjx4664oor1K5duwo/Tt++ffXQQw9p/Pjx+u6771RcXKwmTZroggsu0KBBg0L/iUUgk5Tsb7pWLZV05D4uGi4AAADgKDHTcA0cOFADBw6skte1a9dO999/f6U/Viwx7bvKBhqupdIl1zhOBAAAAESeuFwWHqfPdOgSLNavls077C4MAAAAEKFouHBqzmwppflv0JXnSWtWuM0DAAAARCAaLpwS40vwb4J8hF21xF0YAAAAIELRcOGUmfZdA8el93MBAAAACKLhwikz7cvcx7Vtk+z+Pe7CAAAAABGIhgunrl5DqU79QGm/zXIYBgAAAIg8NFw4ZcaY8me5VnJZIQAAAFAWDRdOT5mGy65aKmutwzAAAABAZKHhwmkx7ToHi/17pJ1b3YUBAAAAIgwNF06LqVFLatoyULNaIQAAABBEw4XTZtoHz3JZ7uMCAAAAAmi4cNrK7sel1ctkS0qcZQEAAAAiCQ0XTt9ZHaSERP9xXq60aZ3bPAAAAECEoOHCaTMpqVKrzEBtVy5xFwYAAACIIDRcCImylxWycAYAAADgR8OFkCi3AfK6VbIFBe7CAAAAABGChguh0eIsKSXVf1xcLK1b6TYPAAAAEAFouBASJjFRyuwYqFkeHgAAAKDhQgiVvayQ+7gAAABwOmxxkUpefFR27SrXUU4LDRdCptx+XJvXyx7KcZYFAAAA0c1+8G9p0Wx5//s7eRPekbXWdaRTQsOF0Gl8plTzDP+xtdLqZW7zAAAAICrZlYtlP/vAX3ietH+vjDFuQ50iGi6EjDFGpl3nQM19XAAAAKgse/CAvL8/HRxodKbMmJ85y3O6aLgQWmX34/qWhgsAAAAnz1or7x/PSAf2+QcSk+S75R6ZatXcBjsNNFwIKdM+eIZLu7bLZu90FwYAAABRxX4+SVq2IFCbK38i07Slw0Snj4YLIWVq15MaNAnUrFYIAACAk2E3rpV9/+/Bgc69ZAaPdBcoRGi4EHJll4cXDRcAAABOwOYekvfS41JxsX+g5hny3fSrqF0ooywaLoRcuf24vs2S9TyHaQAAABDJ/PdtPS2V3opifPL9/B6ZjJpOc4UKDRdCr20nyRz51jp4QNq20W0eAAAARCz72QfS0nmB2lx2rUzbTg4ThRYNF0LOpKVLzVsHapaHBwAAQEXs6uWy498MDnTqKXPhFe4ChQENF8Lih5cVAgAAAGXZ7J3yXnrMv7GxJNWuJ99P75bxxVaLElufDSJGuYUz1iyXLS5yFwYAAAARxRbky3v+L9KhHP9AYpJ8t/5WJr2G22BhQMOF8GjTXkpK9h8X5Evr17jNAwAAgIgQWCRjy/eBMXP9HTItM51lCicaLoSFSUr2N11HsB8XAAAAJMlOelda+E2gNsMula//YIeJwouGC2Fj2ncNHNtvabgAAADinV08R/ajt4MDHbrJXHGTszxVgYYLYWPadw4W61fL5h12FwYAAABO2a0b5b3+VHCgfmP5bvmNTEKCu1BVgIYL4dOslVQ93X/sedKaFW7zAAAAwAl7KEfecw9LBXn+gZRU+X75e/92QjGOhgthY3wJUrvgWS67aom7MAAAAHDCFhfLe/mvUvZO/4Ax8v38HplGZ7oNVkVouBBWZS8rZD8uAACA+GKtlf33C1KZnwPN5dfLdO7lMFXVouFCWJVdOENbN8oe2OcsCwAAAKqWnfSu7NfTArXpfb7MhVc4TFT1aLgQXvUbSbXrBUqWhwcAAIgP3jfTy69I2KaDzE2/kjHGXSgHaLgQVsaY8qsVsjw8AABAzLOrlsq++VxwoEET+e64379Xa5yh4UL4ld2Pa9VSWWvdZQEAAEBY2a0b5b34qFRS4h/IqCnfXQ/IpNdwG8wRGi6EXbkzXHuzpZ3b3IUBAABA2Nh9e+Q985BUuv9qcrJ8d/4/mXoN3QZziIYLYWdqnCE1aR6ouY8LAAAg9tj8w/L+9idpX7Z/wPjk+/lvZFpmug3mGA0XqoRp3yVwbLmPCwAAIKbY4mJ5Lz0ubdkQGDM//rlM1z4OU0UGGi5UibINl77NkvVK3IUBAABAyFhrZd96UVqxODBmhl8u36CRDlNFDhouVI3Ms6WEBP/x4Vxp43q3eQAAABASdtI42a+mBmrT4xyZK250mCiy0HChSpiU6lKZ63e5rBAAACD6ed98LvvRW8GBNu1lfvZrGR9tRilmAlWm3H1cLJwBAAAQ1fx7bT0bHKjfWL47fh+Xe20dDw0Xqowpsx+XvlspW1jgLAsAAABOnX+vrcfYa+sk0HCh6rTMlKql+o+Li6S1q9zmAQAAQKXZ/Xvk/e0hKS/XP5CcLN8v/yBTv5HbYBGKhgtVxiQm+hfPOIL7uAAAAKJLYK+tvaV7bRn5fn6PTKu2boNFMBouVKly93GtpOECAACIFrakRN7Lf5U2l9lr60c/l+na12GqyEfDhSpVbj+uTetkcw+6CwMAAICTYq2V/c/L0vJFgTEz/DL5Bo9ymCo60HChajVpLmXU9B9bK61e5jYPAAAATshO/Uh21pTgQPf+Mlfc5CxPNKHhQpUyxrA8PAAAQBSxi+fIvv+P4EDLTPnYa+ukMUuoeuXu41riLgcAAACOy37/nbzXnvBfmSRJderL98vfyyRXcxssitBwocqV249r13bZ3TucZQEAAEDF7J5d8p57WCos9A+kpsn3qz/K1DjDbbAoQ8OFKmfq1JMaNgnUdtUSd2EAAABwFHs4V96zf5YO7PMPJCTId+tvZRo3cxssCtFwwYmyZ7m4rBAAACBy2OJi//LvWzcGxsy1t8l06OouVBSj4YIT5uxuwWJVlqxX4i4MAAAAAuy416SViwO1uegK+c4b7jBRdKPhghttO0oJCf7jw4ekjevc5gEAAIC8Lz+TnTE5UJse58hcdr3DRNGPhgtOmJTqUqu2gdquWHycZwMAACDc7LpvZd9+KTjQvI3MT+9m+ffTxOzBGdMheFkhC2cAAAC4Y/fvlffiY1JxsX8go6Z8t/+O5d9DgIYLzpS78XLdt7L5h51lAQAAiFe2qEjeS49JB/b6BxIS5Lv1Ppna9dwGixE0XHCnRRupepr/uKREWr3cbR4AAIA4Y62V/c/L0rpvA2PmRz+XyTzbYarYQsMFZ4wvQWrfJVCzPDwAAEDVsrOmyH75WaA25w2XGXCRw0Sxh4YLTpW9rNCuZOEMAACAqmLXrJB955XgQKu2Mj/+hYwx7kLFIBouOFV2A2Tt2Cq7Z7ezLAAAAPHC7s3237dVcmQv1Jq15bvtPpmkJLfBYhANF5wy9RpK9RsFas5yAQAAhJctLpb36v9KBw/4BxIT/c1WrTpug8UoGi44V3Z5eK1a6i4IAABAHLAf/ltauypQmx//QqZ1O4eJYhsNF5wrfx/XElmvxF0YAACAGGaz5st+Oj5Qm36DZM4b7jBR7Et0HcC1NWvW6OOPP9bq1at16NAhpaSkqGXLlho+fLj69u1b4WtmzZqlKVOmaMuWLUpMTFRmZqZGjx6ttm3bVnH6GNG2k+TzSZ4n5R6UNq2XWpzlOhUAAEBMsXt2y/v708GBRmfKXHsbi2SEWVw3XLNnz9bTTz8ta61at26ts88+W/v27dOKFSu0fPlyXXrppbr22mvLveaNN97QpEmTlJycrM6dO6uoqEhZWVlaunSpxo4dq969ezv6bKKXqZ4mtWobOLVtVy6RoeECAAAIGVtSIu+1J/y/3Jak5Gry3fpbmWopboPFgbhtuEpKSvT666/LWqu7775b/fv3Dzy2Zs0aPfTQQ/r44481ZMgQNWzYUJK0fPlyTZo0SRkZGXr44YfVqFGjwPMffPBBvfDCC+rQoYPS09OdfE7RzLTvKlum4dKIq9wGAgAAiCF20rjy921de5tM42YOE8WPuL2Ha+vWrcrJyVGTJk3KNVuSlJmZqS5dushaq/Xr1wfGJ0yYIEkaPXp0oNkqff6wYcN0+PBhzZgxo2o+gRhjzi6zcMbaVbIF+e7CAAAAxBC7dpXsxHcDtek7SL7+gx0mii9x23AlneQeA6VnqwoLC7V8+XJJqvDertKxhQsXhihhnGlxlpSa5j8uKZbWLHebBwAAIAbYw7nyXvs/yXr+gXoNZa75hdtQcSZuG64GDRqoQYMG2rp1q7755ptyj61Zs0ZLly5V/fr11aFDB0nStm3bVFRUpBo1aqhOnaP3KGjZsqUkaePGjeEPH4NMQoLUrlOgtiuXuAsDAAAQI+zbL0l7dvkLn0++m/9HJrW621BxJm7v4fL5fLr99tv1+OOP6+mnn9aECRPUoEED7du3T99++63atGmjO++8U4mJ/inKzs6WpAqbLUlKSUlRWlqacnNzlZeXp9TU1Cr7XGKF6dBNdvEcSZJdwQbIAAAAp8Ob/6Xs3FmB2lz8Y5lWrKpd1eK24ZKk9u3b68EHH9QTTzyhdevWad26dZKk1NRUderUSWeccUbgufn5/nuKkpOTj/l+1apVU25urvLz80+q4Ro7duxRY8nJyXrsscckSXXr1q3U5xNKpY1mvXr1quxjFp87WNlvvegvtm9WbWOVULd+lX38cHMxp7GOOQ095jT0mNPQY05DjzkNPddzWrJvj7L/83KgTurQVbVvuM1/VVGUcj2npyquG66vvvpKL774os466yzdfffdatq0qfbt26cJEyZo/PjxWr58uR588EElJibKWitJx92noPQ5ODWJjZoqoUETlezcKkkqWDpf1YeMdJwKAAAgulhrlfPCY7IHcyRJJiVVNX/1h6hutqJZ3DZc27dv1/PPP69atWrpvvvuU0qKfw+CRo0a6ZZbbtG+ffu0cOFCzZw5U0OHDg2csSooKDjmexYWFkpS4L1O5Mknnzzu49nZ2c6auNLfHOzevbtKP67XtpN0pOE6OGeWcjvHzr5mruY0ljGnocechh5zGnrMaegxp6Hnck69b6bLzv8qOHDlT7QvIVmK8q+vqzk1xpRbobyy4nbRjK+//lolJSXq0qVLhQ1Sv379JEkrVqyQFLy8b8+ePRW+X35+vnJzc5WWlsb9W6eh7PLwdtUSWa/EYRoAAIDoYvfuln3n1eBAh24y51/gLhDit+Hau3evJKl69YpXaSkdP3TokCSpcePGSkpKUk5OToVN14YNGyRJzZqxgdxpad9FKj3dfeig9P1at3kAAACihLVW3r9ekPIO+wdSq8t34y+Pe0sMwi9uG65atWpJUmChjB9au9b/g37pqcvk5GR17NhRkjRnzpyjnl861qNHj1BHjSsmtbrUun2gtssXOUwDAAAQPey8L6TlwT1hzdU/l6kdXQtMxKK4bbh69uwpSVq1apU+++yzco+tWbNGkyZNklR+k+ORI/0LOIwfP17bt28v9/xp06YpNTVVgweza/fpMh27B47tChouAACAE7EHc8pfSnh2N5n+/FwaCeJ20YxWrVrp4osv1oQJE/Taa6/p008/VZMmTbRv3z6tWbNG1loNHTpUnTt3Drymc+fOGjFihCZPnqx7771XnTp1UklJibKysuR5nu68806lp6c7/Kxigzm7u+z4N/3Fhu9kD+XIpNdwGwoAACCC2fdelw75VyVUcjX5rrudSwkjRNw2XJJ0/fXXq23btpo6darWr1+vbdu2KSUlRR06dNCQIUN07rnnHvWam266SS1atNCUKVO0bNkyJSQkqGPHjrriiivUrl07B59FDDqzpVTzDOnAPsl6siuXyPQ+33UqAACAiGRXLJadPSNQm8uvk6nbwGEilBXXDZck9e7dW717V27p8YEDB2rgwIHhCQQZY/xnub6Z7h9Yvkii4QIAADiKLSiQ96/ngwMtzpIZPMpdIBwlbu/hQoT7wX1c1vMchgEAAIhMdvJ70p5d/iIhQb4bfinjY4PjSELDhYhk2neRzJFvz5z90pbvXcYBAACIOHbHVtnPxgdqM/RSmTNbOkyEitBwISKZ9BpSy7MCtS2zxCkAAEC8s9bK+88rUnGxf+CMujKjrnYbChWi4ULEMh2De5qxPDwAAEAZi76RVi4OlL6rb5ZJSXUYCMdCw4WIVXY/Lq37VvZwrrswAAAAEcLm58l757XgwNndpO793AXCcdFwIXI1by2lZ/iPS0qkb7Pc5gEAAIgAdvJ70v49/iIxUb4f/4I9tyIYDRcilvElyHToFqi5jwsAAMQ7u3uH7NSPArUZPlqmQWOHiXAiNFyIbD+4j8ta6zAMAACAW/a/b0jFRf6iVm2ZEVe6DYQTouFCRDNndw0We7Ol7ZudZQEAAHDJrlkuu/DrQG1G3yhTLcVhIpwMGi5ENFPjDKlZ60Btl7NaIQAAiD/WK5H3bpmFMlpmyvQZ4C4QThoNFyJe2dUKuY8LAADEIzt7hrRpfaD2XX2zjI8f5aMBXyVEvLL7cem7FbIF+e7CAAAAVDFbkC/7wb8Dtel9vkzrdg4ToTJouBD5WrWVUtP8x8XF0uplbvMAAABUITvtY+nAXn+RlCwz+ka3gVApNFyIeCYhQWrfJVBzHxcAAIgX9mCO7KfjA7UZPEqmTj2HiVBZNFyICtzHBQAA4pGdPE7KO+wvqqfLXMQy8NGGhgtRwZwdbLi0e4fsrm3uwgAAAFQBu3uH7IzJgdqMuEomLd1hIpwKGi5EBVO7rtSkeaC2yzjLBQAAYpv96C2ppNhf1K4rM3ik20A4JTRciBplVyu0WQscJgEAAAgvu3mD7NxZgdpceq1MUrLDRDhVNFyIGqZzr2CxZpls/mF3YQAAAMLI++itYNGkuUzfgc6y4PTQcCF6tG4nVT9y3XJxsbRyqds8AAAAYWA3fCctnReofZdeK+NLcJgIp4OGC1HDJCT84LLCecd5NgAAQHTyPi5zdqt5G6lrH3dhcNpouBBdOvcMHNqsBbKe5zAMAABAaNm1q6Qye476Lr1GxhiHiXC6aLgQVUzHHpLvyLftwQPSxrVuAwEAAIRQuXu3WrWVylzdg+hEw4WoYtLSpTYdArXNmu8wDQAAQOjY1cukb7MCte/Sazm7FQNouBB1yq5WSMMFAABihTfhnWCRebbUvou7MAgZGi5EnXLLw29aL7s3210YAACAELDfrZRWLwvUvks4uxUraLgQfRo2keo1DJR2GZsgAwCA6OZNejdYnNVBpm1Hd2EQUjRciDrGGC4rBAAAMcNu+E5asThQ+0Zd7TANQo2GC1HJdOkdLFYtlS0ocBcGAADgNHiTxwWLlplS+67OsiD0aLgQnc7qIKWk+o+LCqXVWcd/PgAAQASyWzZIS+YGat/IMdy7FWNouBCVTGKSdHa3QG2XclkhAACIPnbSe8GiaUup7OJgiAk0XIhaP7yPy1rrMA0AAEDl2B1bZRd+Hah9ozi7FYtouBC1TMceUulfSvv3SJs3uA0EAABQCXbqh1LpL4wbNpG69XOaB+FBw4WoZWrUklq1DdSsVggAAKKFPbBP9pvPA7W5YLSMjx/NYxFfVUQ106ln4JiGCwAARAs7fYJUXOQvataW6TPQaR6EDw0XoprpUubG0u+/k83Z5y4MAADASbD5h2VnfhKozdCLZZKSHCZCONFwIbo1aSHVrus/tlZ22UKncQAAAE7EfvGZlJfrL1Kry5x/odtACCsaLkQ1Y8xRqxUCAABEKltcJDv1o0Btzr9Qpnqaw0QINxouRD3TuXewWLFEtqjIXRgAAIDjsPO+9K+uLEmJiTJDL3YbCGFHw4Xo166TlFzNf1yQJ63OcpsHAACgAtba8me3+gyUqVXHYSJUBRouRD2TlCx16Bao7ZK5DtMAAAAcw+pl0pbgvqFm2GXusqDK0HAhJphufQLHdsk8Wc9zmAYAAOBo3rSPg0WHbjJNmrkLgypDw4WYYDr3ksyRb+cDe6Xvv3MbCAAAoAy7c5tUZnEv37BLHKZBVaLhQkww6TWkzLMDtV0yx2EaAACA8uz0jyVr/UWjM6Wzu7sNhCpDw4WYYbqWuaxwMfdxAQCAyGBzD8l+PT1Qm6EXyxjjMBGqEg0XYkbZhks7tshu3+IuDAAAwBH2y0+lwgJ/kZ4h03eQ20CoUjRciBmmbgPpzJaBmssKAQCAa7akRHbGpEBtzr9IpnQ7G8QFGi7EFNO1b+DYLqbhAgAAji2dK+3N9h8nJMgMushtHlQ5Gi7EFNMt2HBpwxrZ0p3cAQAAHPA+L3N2q3t/NjqOQzRciC1NW0h16gdKNkEGAACuFG1a79/s+AgzeJTDNHCFhgsxxRhT7iwXlxUCAABXDk/+b7Bo1kpq3c5dGDhDw4WYU+6ywm+zZA/luAsDAADikpd7SPkzPwnUZtBIloKPUzRciD1t2ks1avmPPU926TyncQAAQPzJ+3ySbH6ev0jLkOl9vttAcIaGCzHH+BLKX1a48BuHaQAAQLyxnqfDnwQvJzTnDmMp+DhGw4WYZLr3DxarlsjmHXYXBgAAxJdVS1WybbP/2BiZgSwFH89ouBCbMjtK6Rn+4+Ji2az5bvMAAIC44c0K3rulzr1k6jZwFwbO0XAhJpnERJkufQK1XcRlhQAAIPzsvj1SmfvHfQM4uxXvaLgQs0yPMpcVLl8oW5DvLgwAAIgL9qupkudJkhLqN5LO7uo2EJyj4ULsatdFSk3zHxcWSssXus0DAABimi0pkf3ys0CdOvxSGV+Cw0SIBDRciFkmKUmmS69AzWqFAAAgrJYtkPZl+48TEpQ6ZJTbPIgINFyIaWVXK7RZC2QLCxymAQAAscybNSVwnNJngBLOqOMwDSIFDRdi29ndpGqp/uOCPGn5Ird5AABATLK7d0grgj9npF54ucM0iCQ0XIhpJrmaTJfegdou+MphGgAAEKvsl59J1vqLBk2U3KmH20CIGDRciHmm17mBY7t0HqsVAgCAkLIlJbLffB6ozfnDZYxxmAiRhIYLse/s7mVWKyyQzVrgNg8AAIgtyxdJB/b6jxMSZfoNdpsHEYWGCzHPJCXJdC2zCfKCLx2mAQAAscb7amqw6NJbJqOmuzCIODRciAum13nBYtlC2fzD7sIAAICYYXP2ScvmB2rfucMcpkEkouFCfGjfRUrL8B8XFcoumec2DwAAiAl29kyppMRf1Kojnd3VZRxEIBouxAWTmCjTvV+gZrVCAABwuqy1sl9PC9Sm/xAZX4LDRIhENFyIG6ZncLVCLV8ke/iQuzAAACD6rV8tbd8cKM05QxyGQaSi4UL8aNtJKr2JtaRYdvEct3kAAEBUs2UXy2jbSaZ+I3dhELFouBA3TEKCTI/+gdrOneUwDQAAiGY2P092fvAWBXPuUIdpEMkSXQeIBPv379eHH36oxYsXKzs7W8nJyapfv746deqk66677qjnz5o1S1OmTNGWLVuUmJiozMxMjR49Wm3btnWQHpVh+gyQnfmJv/g2S3b/HpladdyGAgAAUccu/FoqyPMXqWky3fsf/wWIW3F/hmvNmjX69a9/rcmTJyshIUE9e/ZUZmamDh06pIkTJx71/DfeeEPPP/+8Nm/erE6dOqlNmzbKysrSAw88oHnzWPku4rVuL9Wp7z+2VnYee3IBAIDKs1+VWSyj93kyydUcpkEki+szXHv37tWjjz6qoqIi3XPPPerdu3e5x9euXVuuXr58uSZNmqSMjAw9/PDDatTIf53umjVr9OCDD+qFF15Qhw4dlJ6eXmWfAyrHGOM/yzX5PUlHLiscfpnbUAAAIKrYHVuktSsDtWHvLRxHXJ/hevvtt5Wbm6vrrrvuqGZLktq0aVOunjBhgiRp9OjRgWZLkjIzMzVs2DAdPnxYM2bMCG9onDbTZ0Cw2LROdvsWd2EAAEDUKXt2S01bSM3bHPO5QNw2XIcOHdLs2bNVvXp1DR48+ITPLyws1PLlyyVJffv2Perx0rGFCxeGNihCzjRuJjVrFajt3JnuwgAAgKhii4tlZ38eqM05Q2WMcZgIkS5uLylcvXq1ioqK1KlTJyUmJmrOnDn69ttvVVxcrCZNmqhfv36qVatW4Pnbtm1TUVGRatSooTp1jl5koWXLlpKkjRs3VtWngNNg+gyQ3bRekv+yQnvptfxlCQAATmz5Qilnv/84MVGm70CXaRAF4rbh2rzZv0ldzZo19cc//lFr1qwp9/jbb7+t22+/Xf369ZMkZWdnS1KFzZYkpaSkKC0tTbm5ucrLy1NqamoY0+N0mV7ny77/T8laKXuntO5bqU1717EAAECE874us1hG174y6TUcpkE0iNuGKzc3V5L0xRdfKDExUbfeeqt69uyp/Px8TZkyRRMnTtSzzz6rxo0bq3nz5srPz5ckJScnH/M9q1WrptzcXOXn559UwzV27NijxpKTk/XYY49JkurWrXsqn1pIJCb6vzXq1avnLENY1aunvR27q3CZ/xLQlKy5qtHv/LB+yJifUweY09BjTkOPOQ095jT0mNOTU7Jvj3ZnLQjUtUZeoWrHmDPmNPSidU7j9h4uz/MkSSUlJbrxxhs1ePBg1ahRQ/Xr19cNN9ygvn37qri4WB999JEkyVorSce97Kz0OYgOKQMvDBznfTVNtqjQYRoAABDp8mdNkbwSSZKvbgMld+7lOBGiQdye4So9A2WM0YABA456fNCgQZozZ45WrlxZ7vkFBQXHfM/CQv8P7CkpKSeV4cknnzzu49nZ2c6auNLfHOzevdvJx68K9qxOUnKyVFgoezBHu6d/ItMjfJsWxsOcVjXmNPSY09BjTkOPOQ095vTklEydEDi2fQYoe+/eYz6XOQ09V3NqjCm3Qnllxe0ZrtIvWK1atZSUlHTMxw8cOCApeHnfnj17Kny//Px85ebmKi0tjfu3ooRJrV5uV3jvm+kO0wAAgEhmN2+QtgYXRzP9TrzKNSDFccNVuqpgbm5uhWeRDh06JCl4tqpx48ZKSkpSTk5OhU3Xhg0bJEnNmjULV2SEgek/JFgsXyh7YJ+7MAAAIGLZOWX2Wm2ZKdOgsbswiCpx23A1a9ZM9evXV2Fhob777rujHl+xYoUkqVUr/35NycnJ6tixoyRpzpw5Rz2/dKxHjx7hioxwaNtJqn1kcRLPk50z02kcAAAQeaxXIjv3i0Bt+g1ymAbRJm4bLkm69NJLJUn/+Mc/lJOTExhfv369Jk6cKEkaNmxYYHzkyJGSpPHjx2v79u2B8TVr1mjatGlKTU09qU2UETmMz1fukgD7zXQWPwEAAOWtypIOHLlfKyFBpud5bvMgqsTtohmSNGTIEC1btkxz5szR3XffrczMTBUUFGj16tUqLi7WkCFD1Ldv38DzO3furBEjRmjy5Mm699571alTJ5WUlCgrK0ue5+nOO+9Uenq6w88Ip8L0HyI7aZy/2LZJ2rhWanGW21AAACBilLucsGMPmQz23sLJi+uGy+fz6e6779bUqVP1+eefBy4jbN26tYYNG6bzzz96X6abbrpJLVq00JQpU7Rs2TIlJCSoY8eOuuKKK9SuXbuq/hQQAqZ+I+msDtJ3/hUp7dfTZWi4AACAJJufJ7todqD2cTkhKimuGy7J33RdcMEFuuCCC076NQMHDtTAgQPDFwpVzvQfIlvacM2bJTvmpzJJx97kGgAAxAe7eI5UeGRboNQ0ib23UElxfQ8XUMr0PEdKruYvDuf6/3IFAABxr+zlhKbnOfxCFpVGwwVIMinVZXqeG6jtF586TAMAACKB3b/Hv2DGEaYvlxOi8mi4gCPM+WUuK129THbHVndhAACAc3buF5L1/EWd+lKb9m4DISrRcAGlWrWVmjQPlParzxyGAQAArpW7nLDvQBkfPzqj8viuAY4wxsicFzzLZb+eLltU5DARAABwxW7ZIG35PlBzOSFOFQ0XUIbpO1AqvRn2UI7skrlO8wAAADfs7JnBomWmTMMmzrIgutFwAWWYtHSZHucEavsli2cAABBvrFciO29WoDbsvYXTQMMF/EC5xTNWLZXdtd1dGAAAUPW+zZL27/UfJyTI9DzPbR5ENRou4IfatJcanRkoWSIeAID4YmcHF8tQxx4yGTXchUHUo+ECfsAYI3P+8EBtv54qW7rDPAAAiGk2P0920exA7eNyQpwmGi6gAqb/ECm5mr84dFB2/lduAwEAgCphl8yVSn/Rmpomde7lNhCiHg0XUAFTPb3c8q/284my1jpMBAAAqoKd90Xg2PToL1O6ejFwimi4gGMwg0YEi03rpPWr3YUBAABhZw/lSCsXB2rT+3yHaRAraLiAYzBNW0htOwVq+/lEd2EAAEDY2UXfSCUl/qJmbaltR7eBEBNouIDj8A0aGTi2C7+WLV0iFgAAxBw7t8zlhD3PkfElOEyDWEHDBRxP1z7SGXX9xyUlsl9+5jYPAAAIC7tvj/TdikDN5YQIFRou4DhMQoLMgAsDtZ31iWxxkcNEAAAgHOz8L6XSBbLqNZRaZroNhJhBwwWcgDn/AikxyV8c2Fdu9SIAABAbyq1O2Os8GWMcpkEsoeECTsBk1JQps+mh/exDlogHACCG2J3bpI1rA7XpPcBhGsQaGi7gJJhhlwaLrRulFYvchQEAACFV7uqVJs1lmjRzFwYxh4YLOAmm0Znldpr3PvvQXRgAABAy1trylxOyWAZCjIYLOEm+4ZcHi1VLZTetcxcGAACExuYN0o4tgdL0Os9hGMQiGi7gZGWeLbU4K1BaznIBABD1yl1O2KqtTL2G7sIgJtFwASfJGCNT5iyXnf+l7J7dDhMBAIDTYT3Pvxz8ESyWgXCg4QIqwXTvJ9Wp7y88T/azD9wGAgAAp27dt9LeI788NT6ZXue4zYOYRMMFVIJJSJAZflmgtl9+Jntgn7tAAADglJW7nLBdJ5kaZ7gLg5hFwwVUkjl3mFTzyF/IRYWc5QIAIArZkhLZhV8HalYnRLjQcAGVZJKrlb+Xa+YnsgcPOEwEAAAqbdVSqfTf78RE/20DQBjQcAGnwAy4UEqv4S8KC2SnfuQ2EAAAqJRylxN27CFTPd1dGMQ0Gi7gFJhqKeXv5ZoxSTb3oLtAAADgpNmiQtnFswM1qxMinGi4gFNkBo2QSn8blp8nO22C20AAAODkLF8k5ef5j6ulyHTu5TYPYhoNF3CKTEp1maGXBGo7/WPZgzkOEwEAgJNRbu+tLr1lqlVzmAaxjoYLOA1myKjgWa68w7JT3ncbCAAAHJctKJDNmh+oTc9zHaZBPKDhAk6DqZ4uc9EVgdp+Pkl2b7bDRAAA4LiWL5AK8v3HKalSx+5u8yDm0XABp8kMGiXVrO0viotkJ77jNhAAADgmO/+rwLHp2kcmKdlhGsQDGi7gNJlq1WRGXR2o7dfTZHdsdZgIAABUxBbkyy7jckJUrZA3XHfeeac+/PBDHTjARrCIH+bcYVK9hv7C82Q/esttIAAAcBSbNV8qLPQXqWlSh25uAyEuhLzh2rVrl/7zn//otttu05NPPqmsrKxQfwgg4pjERJlLrw3UdsFXshvWOEwEAAB+yC744eWESQ7TIF6EvOG6/PLLVbt2bZWUlGju3Ln6y1/+ojvvvFMfffQRZ70Q00yv86QzWwZq751XZa11mAgAAJSy+YelZQsDtel1nsM0iCchb7h+9KMf6fnnn9e9996rbt26yefzadeuXXr77bd122236amnnuKsF2KS8fnkG/Oz4MD61bJzZ7kLBAAAAuzS+VLRkcsJq6dL7Tu7DYS4kRiON/X5fOrRo4d69OihvXv36vPPP9eMGTOUnZ2tOXPmaM6cOapfv76GDh2qgQMHqmbNmuGIAVQ5066z1L2ftGi2JMn+9w3Zbn1lqqU4TgYAQHwrt9lx934yiVxOiKoR9lUKa9eurSuvvFLPPfecfve736lXr15KSEjgrBdilu/Kn0iJR36XsX+P7JT/ug0EAECcs4dzpRWLAjWrE6IqheUMV0WMMeratau6du2qvXv36m9/+5tWrVqlkpKSwFmvhg0b6uKLL9bgwYPl87FiPaKTqddQZthlsp+8L0myn34ge+4wqV49x8kAAIhPdslcqbjYX6RnSO24nBBVp0q7muzsbI0bN06///3vtWrVqsB4ixYt5PP5tGPHDr366qv6/e9/r5ycnKqMBoSUGXGlVPMMf1FUKO+dV90GAgAgjpVbnbB7f5mEBIdpEG/CfobL8zwtWLBA06dPV1ZWljzPkySlp6dr4MCBGjZsmBo2bKj9+/frs88+06RJk7R+/Xq9/fbbuvXWW8MdDwgLk1JdZvSNsv942j+wZK7y58xUSt+BLmMBABB3bO4haeWSQM3lhKhqYWu4du3apenTp2vmzJnav39/YDwzM1PDhw9Xv379lJgY/PC1atXSmDFj1KNHD91///1avHhxuKIBVcL0GyT7zXRp9TJJUs4r/6fkTj0dpwIAIL7YJXOkkiOXE2bUlDI7ug2EuBPyhmvOnDmaNm2ali9fHtiDKDU1Veeee66GDx+uZs2aHff1rVu3Vq1atco1aUA0MsbId93t8h76lVRcJG9vtg699ZJ0+Y2uowEAEDfKrU7Y4xwuJ0SVC3nD9dRTTwWOW7RooWHDhuncc89VSsrJL4td9swXEM1MwyYyI6+S/ehtSdLhT8bL17mPTOt2jpMBABD77KEcadXSQM3lhHAh5J1NUlKS+vfvr+HDh6tNmzan9B7PP/98iFMB7pgLr5Cd96W0fbNkrbw3n5PvD0/KJCW7jgYAQEyzi2ZLR9YPUM0zpLPauw2EuBTyVQpffvll3X777afcbAGxxiQmyXf9HcGBbZtkP/y3u0AAAMSJcqsT9jhHxsflhKh6IW+40tLSQv2WQNQzZ3VQ9VFjArWd+pHst2z2DQBAuNiDB6RvlwVqLieEK+wuDFSRjOtvU0LTFv7CWnn/eFr28CGnmQAAiFV24TeSPXI5Ya06EvdPwxEaLqCKmGopqvXrB6XS1ZH2Zsv+5xWnmQAAiFXlLifseY6Mjx974QbfeUAVSmrdVubiHwdqO2emvNkzHCYCACD22AP7pDUrAjWXE8IlGi6gipkLryh3WYP99/OyW753FwgAgBhjF5W5nLB2PalVW7eBENdouIAqZhIS5Pv5PVJahn+gsFDei4/JHs51GwwAgBhRbrPjnufKGOMwDeIdDRfggKlTX76bx0ql/wDs2ibvn8/IWus2GAAAUc7u2yOtXRWouZwQrtFwAY6Yjj1kRv0oOLB4juwn77sLBABADLALv5ZKf4FZt4HUgr1h4RYNF+CQGXW11LF7oLYf/KvcqkoAAKByyq9OyOWEcI+GC3DI+Hzy/WysVK9hYMx7/SnZMpdCAACAk2P37pbWfRuoTS8uJ4R7NFyAYya9hny/+qNUPd0/UFwk7/mHZXdtcxsMAIAoYxd8HSzqN5LObOUuDHAEDRcQAUzDpvLdcb+UmOgfOHRQ3jMPye7f6zYYAABRhMsJEYlouIAIYTI7ytx0V3Bg13Z5T/4/2Zz9zjIBABAtbPZOacOaQM3lhIgUNFxABPH1GSAz+sbgwPbN/qbrUI67UAAARIFyi041bCI1aeEsC1AWDRcQYXwXXSFzcZnl4rdulPfUH2VzD7oLBQBAhCt7/5bpeR6XEyJi0HABEchc/GOZi64MDmxaL+/x+/yrLwEAgHLsrm3SxrWBms2OEUlouIAIZIyRufx6meGXBwe3b5b32G9lt25yFwwAgAhUbnXCxs1kmjRzFwb4ARouIEIZY2SuvElmVJnLC/dly/vrfbLfrXQXDACACPPD1QmBSELDBUQwY4x8l14jc+1tUum16IcPyfu/P8ib+YmstW4DAgDgmN2xVdq8IVDTcCHS0HABUcA38CL5bv2tlJjkHygpln3rRdk3/iZbVOg2HAAADpVbnbBpC5lGTd2FASpAwwVECdO9v3z/87BUs3ZgzH493b+Yxo4tDpMBAOCOnf9l4JizW4hENFxAFDFt2sv3hyelNu2DgxvXyvvz3fI+nyjree7CAQBQxezWTdK24GJSbHaMSJToOkAkOXTokO6++27l5OSocePGevrpp4/53FmzZmnKlCnasmWLEhMTlZmZqdGjR6tt27ZVFxhxydSqLd//PCw77nXZGZP9g4WFsv95RXbJXPluvFOmTn23IQEAqAJ2QfDslpq3kanf2F0Y4Bg4w1XGG2+8oYMHT7y57BtvvKHnn39emzdvVqdOndSmTRtlZWXpgQce0Lx586ogKeKdSUyS75pb5bv9fimjZvCBVUvl/fF2eRPf5d4uAEBMs9bKzi+zOiFntxChOMN1xLJlyzRr1iwNHTpU06ZNO+bzli9frkmTJikjI0MPP/ywGjVqJElas2aNHnzwQb3wwgvq0KGD0tPTqyo64pjp1le+1u3k/et5aclc/2BhoexHb8nO/ly+MT+TOveSKV3hEACAWLF5g7Rza6Dk/i1EKs5wSSosLNSrr76qpk2b6uKLLz7ucydMmCBJGj16dKDZkqTMzEwNGzZMhw8f1owZM8KaFyjL1Kgl3+33y/zs1+XPdu3aLu+5h+U9/lvZVUvdBQQAIAzKXU7Yqi2X0yNi0XBJeu+997Rz507dfPPNSkhIOObzCgsLtXz5cklS3759j3q8dGzhwoXhCQocgzFGvr6D5Hv4RZkhF0umzB/tdd/Ke/L/qeT//iC7cjF7dwEAoh6XEyKaxP0lhRs3btTEiRM1cOBAdejQQbt27Trmc7dt26aioiLVqFFDderUOerxli1bBt4TcMFUT5f50c9lzx0q779vSsvLNP/fZsn7Nsu/R8mwy2R6nydTuq8XAADR5PvvpOyd/mNjZHrQcCFyxfUZLs/z9PLLL6t69eq67rrrTvj87OxsSaqw2ZKklJQUpaWlKTc3V3l5eSHNClSGadpSCXc9IN+9j0mZHcs/uOV72X88Le/en8p7/x+yO7ZW/CYAAESosntv6awOMmdU/LMZEAni+gzXlClTtHbtWt1+++3KyMg44fPz8/MlScnJycd8TrVq1ZSbm6v8/HylpqYe9/3Gjh171FhycrIee+wxSVLdunVPmClcEhP93xr16tVzliHWOJnTegNk+52vwuWLdPjDt1Ww8JvgYwcPyH76geynHyjp7K5KPXeYUvoPkq/mGVWX7zTxfRp6zGnoMaehx5yGXjTNqfU87V40W6UXyNcYeJGqR2DuaJrTaBGtcxq3DVd2drbeeecddejQQQMHDjyp15Te+3K8Fd+4PwaRxhijap16qFqnHiratF6HP3pbeV9MlcosG1+0YomKVixRzqtPKrlTD6WcN1QpfQfIl17DYXIAAI5WtHqZvD1HbgHx+VSt/yC3gYATiNuG67XXXlNxcbFuvvnmk35N6RmrgoKCYz6nsND/Q2xKSsoJ3+/JJ5887uPZ2dnOGrjS3xzs3r3bycePRRExp6kZ0o9+Id/F18jOmSn7xafStk3Bx70SFS6dp8Kl85Tz4l+ls7vJ9DhHpksvmbQTnwWuahExpzGGOQ095jT0mNPQi6Y59aZODBZtO2lvUYkUgbmjaU6jhas5NcaUW528suK24Vq0aJHS0tL02muvlRsvKiqS5G92HnzwQUnSfffdp5SUlMAlfnv27KnwPfPz85Wbm6u0tLQTXk4IuGTSMmSGXCw7eJT0/Xey8770L6+7f2/wSSXFUtZ82az5sj6f1LaTTNc+Ml37ytR2d7krACB+Wa9EduHXgZq9txAN4rbhkqTc3FytXLmywscKCwsDj5WUlEiSGjdurKSkJOXk5GjPnj1HLZ6xYcMGSVKzZs3CmBoIHWOM1DJTpmWm7FU/kdaukp3/pf8fs4MHgk/0PGnVUtlVS2X/84rUvI1Mt74y3ftJDZuysTIAoGqsWSEd2Oc/Tkjw/zsERLi4bbjGjRtX4fiuXbv0y1/+Uo0bN9bTTz9d7rHk5GR17NhRixcv1pw5czRy5Mhyj8+ZM0eS1KNHj7BkBsLJ+HxS5tkymWfL/ujn0uplsgu+kl0yt3zzJUkb18puXCv74b+lBk2ONF/9pRZtaL4AAGFTdu8tte8qw73GiAJx23CdqpEjR2rx4sUaP368unfvHriec82aNZo2bZpSU1M1ePBgxymB02MSEqQOXWU6dJW97jZp/WrZxXNlF8+Wdu8o/+SdW2Wn/Fd2yn/9zVffgTJ9BsjUa+gmPAAgJtmSEtlFwdV22ewY0YKGq5I6d+6sESNGaPLkybr33nvVqVMnlZSUKCsrS57n6c4771R6errrmEDIGF+C1KaDTJsOslfeJG3bJLt4tuziudKmdeWfvHOr7EdvyX70lv81fQfK9DxXJo0/EwCA0/RtlnQox3+cmCjTta/bPMBJouE6BTfddJNatGihKVOmaNmyZUpISFDHjh11xRVXqF27dq7jAWFjjJGaNJdp0lwa9SPZPbtkl8z13/P13Q/uh1y7UnbtStl3XvEvtDFopH9zSi45BACcgnKbHZ/dXaZ6mrswQCXQcP1A/fr1j3l/V1kDBw486f27gFhl6tSXGXKxNORi2eydsnNnyc6ZKe3YEnxScbH/XrAFX0lNW8oMHinTd5BMUpKz3ACA6GKLi/yXtR9hep3nMA1QOTRcAELC1G0gM3KM7IirpE3r/Pt8zZ1VfsGNLRtk33xO9uO3ZYZfLnP+BTLVTrxnHQAgzq1YIh3O9R8nJ8t06e00DlAZNFwAQsoY4182vnkb2St/ImXNlzdjkrRqafBJ+/fKjntddvJ7MheOlhk8SiYp2V1oAEBEswvKXE7YqadMCvudInrQcAEIG5OQIHXrq4RufWW3bZKdPlH2m2lScbH/CYdyZN//p+znk2Quv06m9wD/8vQAABxhiwr9W5Qc4eNyQkQZfrIBUCVM42byXX+7fI++KjP8MqnspYR7d8u+/pS8v4yVXb/aWUYAQARatlDKz/MfV0uROvZ0mweoJBouAFXK1Koj31U/le/R12SGXiIllDnRvmm9vMfulffWS7Kl1+oDAOKaXRDc7Nh06SNTrZrDNEDl0XABcMJk1JDv6pvl+/ML5VebslZ25mR5f7xddskcdwEBAM7ZgnzZpfMCNZsdIxrRcAFwytRrKN8tv5HvN49Ijc4MPnBgn7znH5H35nOypZeSAADiis1aIBUW+IvUNOns7m4DAaeAhgtARDCZHeX749Myl10nlVmx0H75mbw/3y27YY3DdAAAF+z8LwLHpmsf9nBEVKLhAhAxTGKSfCPHyPfHp6XmbYIP7Nou7/H75M38RNZaZ/kAAFXH5h32L5hxBJsdI1rRcAGIOKZhU/nu+6vMiDGSOfLXVEmx7Fsvyr7xrGxRoduAAICws0vnSsVF/iItQ2rfxW0g4BTRcAGISCYxUb7Lr5Pvnr9INc8IjNuvp8l7/D7ZfXscpgMAhJudX2Z1wu79ZBLZPhbRiYYLQEQzmWfL94enpNbtgoMb18p77Dcq3rzBXTAAQNjYQznSikWBmssJEc1ouABEPFOrtnz3/EVm4Ijg4N5s7bnvFypcscRVLABAmNgFX0slJf6i5hlS245uAwGngYYLQFQwiUnyXXurzI9vkYyRJNncg9r74F2yC79xnA4AEEp23qzAsel1nowvwWEa4PTQcAGIKr7Bo+S79bdS4pGlgYsK5b3yV3nzvjj+CwEAUcHu2SV9tzJQmz4DHKYBTh8NF4CoY7r3l2/sn2XSM/wDnif72pPy5sx0mgsAcPps2V+gNWhSfpsQIArRcAGISuasDqr98AsyNWr5B6wn+/en5H3zudNcAIDTY+eWuZyw9/kyRy4jB6IVDReAqJXUoo1q//k5KaOmf8Ba2X8+I2/+l26DAQBOid3yvbR1Y6A2fbmcENGPhgtAVEtq3tq/V1fgTJeVff0p2RWLneYCAFRe2bNbapkpU7+xuzBAiNBwAYh6pnEzf9NVek9XSbG8Fx+VXb/abTAAwEmznlfu/i3T+3yHaYDQoeECEBNMozPl+9UDUrUU/0BBvry//Ul2+2a3wQAAJ2ftKmnvbv+x8bHZMWIGDReAmGFaZsp3+/1SQqJ/IPegv+k6eMBtMADACZW7nLB9F5maZ7gLA4QQDReAmGI6dJXv5rGBzZGVvVPeC4/KFhW5DQYAOCZbXCS78OtAbfpwOSFiBw0XgJhjep4rc8WNwYG1K2X/9byste5CAQCObcViKfeg/zgpWaZbP7d5gBCi4QIQk8zwy2XOGRKo7ezPZaeMd5gIAHAs5fbe6tJbJrW6wzRAaNFwAYhJxhiZ626XMs8OjNkP3pRdvshhKgDAD9n8w7JL5wZqLidErKHhAhCzTGKSfLf+TqrX0D9grbzX/092z263wQAAAXbxXKmw0F9UT5c69nAbCAgxGi4AMc1k1JDvtt9JScn+gUMH5b38uGwxi2gAQCSwc2cGjk3Pc2QSk9yFAcKAhgtAzDNntpS57rbgwIY1suP+7i4QAECSZHP2SauWBmrTZ4DDNEB40HABiAu+/kNkzhseqO2MSfLmf+kwEQDAzv9a8jx/Ubuu1KaD20BAGNBwAYgb5se3SM1aB2r7rxdk9+xymAgA4pudMyNwbHqdL+PjR1PEHr6rAcQNk5Qs362/lVJS/QN5ufJef1LWK3EbDADikN26Sfr+u0Bt+g1ymAYIHxouAHHF1Gsoc+2twYHvVspOft9dIACIU3b29GDRvI1Mk+buwgBhRMMFIO74+g6S6R28MdtO+I/s+tUOEwFAfLElJbJzZgZq03+wuzBAmNFwAYhL5tpbpTr1/YXnyXvt/2QLCtyGAoB4sXKxdGCf/zghUaY3mx0jdtFwAYhLpnqafDePlcyRvwZ375D94E23oQAgTthvPg8WXXrJpNdwFwYIMxouAHHLtOkgc8Hlgdp+PlF2zQqHiQAg9tncQ7JL5gRqX/8hDtMA4UfDBSCumUt+LDU6019YK++Nv3FpIQCEkZ3/hVRc7C8yakpnd3cbCAgzGi4Acc0kJct306+Clxbu2i774b/chgKAGFb2ckLTZ6BMYqLDNED40XABiHumVVuZ4ZcFajt9guzaVe4CAUCMsts3SxvWBGpzDqsTIvbRcAGAJHPpNVLDpv7CWnn/fkG29JIXAEBIlFsso1krmaYt3YUBqggNFwDoyKWFN/4yOLB1o+zUj9wFAoAYY70S2TkzArVhsQzECRouADjCtOkgc/4FgdpO/I/s7h0OEwFADFm5RNq/13+ckMDeW4gbNFwAUIYZfaN/1SxJKiyU99aLsta6DQUAMaDc5YSdesmU/l0LxDgaLgAow6Sly1x9c3BgxWLZBV+5CwQAMcAePiS7uMzeWyyWgThCwwUAP2B6ny916Bao7bjXZfPzHCYCgOhm538lFRf5i/QaUscebgMBVYiGCwB+wBgj37W3SolJ/oH9e2UnjXMbCgCimP1meuDY9BkgU/r3KxAHaLgAoAKmfiOZCy4P1HbqR7I7tjpMBADRyW7fIq1fHahZnRDxhoYLAI7BXHSVVLuuvygplvfuayygAQCVZL+eFiyatpRp1spdGMABGi4AOAZTrZp8V/00OLB8oZQ1310gAIgytrio/OWE53B2C/GHhgsAjqfHOVK7zoHSe/c12aIih4EAIIosnS8dPOA/TkyS6TfIbR7AARouADgOY4x8P7pF8h3563L3DtnPJ7oNBQBRwvvi08Cx6d5fJi3DYRrADRouADgB06SZzKCRgdpOGid7MMdhIgCIfDZ7p7RqSaA25w93FwZwiIYLAE6CGXW1VD3NX+Tlyk74j9tAABDh7FdTpdKFhuo3ljI7ug0EOELDBQAnwaTXkBl5daC2sz7xL3UMADiKLSkptzqhOX+4jDEOEwHu0HABwEkyg0ZK9Rr6C8+T9/4/3AYCgEi1fKG0f6//OCFRpt9gt3kAh2i4AOAkmaQk+a64KTiQNV921VJneQAgUnlffhYsuvaWqVHLWRbANRouAKiM7v2kNh0CpfffN2Q9z2EgAIgsds8uKWtBoPadd4HDNIB7NFwAUAnGGPmu+klwYONa2YXfuAsEABHGfvGpZI/8IqpeQ6l9F7eBAMdouACgkkyrtlL3/oHafvgv2eJih4kAIDLY4iLZMpcTmgEXyfj4cRPxjT8BAHAKfJdfF9wMedd22a8+O/4LACAO2EWzpYMH/EViksw5Q9wGAiIADRcAnALTsKnMOUMDtZ3wjmx+nsNEAOCenfVJ4Nj0Ok8mvYbDNEBkoOECgFNkLv6xlJzsL3L2y0772G0gAHDIbt0orVkRqM3Ai9yFASIIDRcAnCJzRh2ZIRcHavvZB7K5Bx0mAgB37Mzg2S01ay21zHQXBoggNFwAcBrMhVdI1dP8Rd5h2c8+chsIAByw+Ydl58wI1GbgRTLGuAsERBAaLgA4DaZ6uszwywO1nT5B9mCOw0QAUPXsN59LpfexpqbJ9B7gNhAQQWi4AOA0mSGjpLQMf1GQJ/vZB24DAUAVsp4nO31ioDbnDJWpVs1hIiCy0HABwGkyKdVlLhgdqO3nE2Vz9jlMBABVaMUiadc2/7ExMoNHus0DRBgaLgAIATNohJRR018UFshOGe82EABUEW/ahGDRpbdMvYbuwgARiIYLAELApKT6F9A4ws78RHb/XoeJACD87PbN0srFgdo39BKHaYDIRMMFACFiBlwk1TzDXxQVyk75r9tAABBmdnqZs1tNW0iZHZ1lASIVDRcAhIipVk3moisDtZ01RXZvtsNEABA+NveQ7OwyS8EPuZil4IEKJLoO4EpBQYGWLl2qhQsXat26ddq9e7c8z1PDhg3Vp08fjRo1SikpKRW+dtasWZoyZYq2bNmixMREZWZmavTo0Wrbtm0VfxYAIo05/wL/ma39e6XiItlP3pe59lbXsQAg5OxXn0mFBf4ivYZMH5aCByoSt2e4vvrqKz3xxBOaMWOGrLXq0qWL2rVrp127dmncuHH63e9+pwMHDhz1ujfeeEPPP/+8Nm/erE6dOqlNmzbKysrSAw88oHnz5jn4TABEEpOULDNiTKC2X34mu2eXw0QAEHq2uEi2zGIZ5vwLZZKSHSYCIlfcnuFKTEzU8OHDNXLkSDVq1Cgwvm/fPj322GPasGGD/vnPf+quu+4KPLZ8+XJNmjRJGRkZevjhhwOvW7NmjR588EG98MIL6tChg9LT06v88wEQOcy5w2SnvC/tzZZKimUnjZO54ZeuYwFAyNh5X0j79/iLxET/Sq0AKhS3Z7gGDBigm2++uVyzJUlnnHGGfvazn0mS5s2bp+Li4sBjEyb4f5MzevTocq/LzMzUsGHDdPjwYc2YMUMA4ptJSpIZWeYs1zfTZXfvcJgIAELHWiv7aXCDd9N3kEyt2g4TAZEtbhuu42nevLkkqaioSAcPHpQkFRYWavny5ZKkvn37HvWa0rGFCxdWUUoAkcz0HyLVqe8vSkpkJ73rNhAAhMryhdK2TYHSDL/cYRgg8tFwVWDnzp2SpISEhMDlgdu2bVNRUZFq1KihOnXqHPWali1bSpI2btxYdUEBRCyTmCRz8Y8CtZ0zk3u5AMQEr+zG7l37yDRq6i4MEAVouCowefJkSVLXrl2VlJQkScrO9i/tXFGzJUkpKSlKS0tTbm6u8vLyqiYogIhm+gwsf5ar7A8pABCF7IY10prlgdp3AWe3gBOJ20UzjmXRokWaMWOGEhISdPXVVwfG8/PzJUnJycdegadatWrKzc1Vfn6+UlNTT/ixxo4de9RYcnKyHnvsMUlS3bp1Kxs/ZBIT/d8a9erVc5Yh1jCnoRcNc3r4qhuV89L/SpLs19NU+4bblFDb3Z/tE4mGOY02zGnoMaehd7Jzuu8fT+nIQvBKatdJdfqxFPyx8H0aetE6p5zhKmPLli169tlnZa3V9ddfrxYtWgQes9ZK0nE39Ct9DgCUSh08Ur4zjjRYRYXK/fg/bgMBwCkq3vK9CmbPDNRpl13rLgwQRTjDdcSePXv0yCOPKDc3V6NGjdKIEeWXNy09Y1VQUFDRyyX5F9aQdMwNk3/oySefPO7j2dnZzpq40t8c7N6928nHj0XMaehFy5zaoZdI7/1dknT4k/HKHzBCJr2G41QVi5Y5jSbMaegxp6F3MnPq/fsVqfTnkoZNldOynQ7yNTgmvk9Dz9WcGmOOWtm8MjjDJSknJ0cPP/ywsrOzNXDgQF1//fVHPaf08r49e/ZU+B75+fnKzc1VWlraSV1OCCB+mAEXSukZ/qIgX3b6hOO/AAAijN21TXberEBtRl4l4+PHSOBkxP2flLy8PD366KPaunWrevfurVtvvbXCywYbN26spKQk5eTkVNh0bdiwQZLUrFmzsGcGEF1MtRSZoZcGavv5RNnDuQ4TAUDl2MnvS57nL+o3kul1vttAQBSJ64arqKhIf/3rX7Vu3Tp16dJFd999t3zH+G1NcnKyOnbsKEmaM2fOUY+XjvXo0SN8gQFELTNohJRa3V8czpWdOdltIAA4STZ7p+ycGYHajLhKJiHBYSIgusRtw+V5np555hmtWLFC7du31z333BNY+eRYRo4cKUkaP368tm/fHhhfs2aNpk2bptTUVA0ePDisuQFEJ1M9XWbQyEBtp30se5x7QgEgUtgp/5VKSvxFnfr+LS8AnLS4XTRjypQpmjdvniQpIyNDr732WoXPu/7661Wjhv/m9s6dO2vEiBGaPHmy7r33XnXq1EklJSXKysqS53m68847AxslA8APmaGXyE77WCoskA4ekP3yU5mhl7iOBQDHZPdmy349LVCbEVfKnOAX1ADKi9s/MYcOHQoclzZeFbnqqqsCDZck3XTTTWrRooWmTJmiZcuWKSEhQR07dtQVV1yhdu3ahTUzgOhmMmrKnH+h7LSPJEn20/GyAy6SObLBOgBEGjt5nFRc7C/OqCvTb4jbQEAUituGa8yYMRozZswpvXbgwIEaOHBgaAMBiAvmgstkZ07y/wCzf6/s7Oky51/oOhYAHMXu2i771dRAbS66kl8QAacgbu/hAgAXTK06MucMDdR2ynjZ0nsjACCC2I/fDt67VbeBzHnD3AYCohQNFwBUMXPBaKl0RdTdO2TnfeE2EAD8gN3yfbm/m8wl18gkcnYLOBU0XABQxUy9hjJ9BgRqO+W/sqX72wBABPA+/Ldkrb9o3EymD/tuAaeKhgsAHDAXXSmVbrK+bZOUNd9tIAA4wq77VloaXFDMd9l1Mj723QJOFQ0XADhgGp0pde0TqL1P3pct/W0yADhirZU3/o3gQMvMcn9XAag8Gi4AcMR30ZXBYv1qac0Kd2EAQJIWzy73d5Hv8utlSs/GAzglNFwA4IhpmSm16xyovSnvO0wDIN7ZokJ57/0jONCpp0z7Lu4CATGChgsAHPJddEWwWL5IdtM6d2EAxLXcCeOk7J3+wueT76qfug0ExAgaLgBwqX1XqXmbQGmnjHeXBUDcKtm/V7llzm6ZQSNlGjV1mAiIHTRcAOCQMabcvVx2wdeyu7Y5TAQgHh166xXZvMP+onq6zMU/chsIiCE0XADgWrc+UoMm/mPrcZYLQJWyG9cqb9rHgdpcco1MWobDREBsoeECAMeML0HmwtGB2s7+XHb/HoeJAMQLW1Ii783ng5scN2wqM+BCt6GAGEPDBQARwPQdKNWq4y+Ki2Wnfnzc5wNAKNgZk6Qyi/X4rr1VJjHRYSIg9tBwAUAEMIlJMsMvC9R21hTZ3EPuAgGIeXbvbtkP3wrUKYMukimzVQWA0KDhAoAIYc4bLpXeN1GQ5//NMwCEifefV6WCPEmSyaihGjfd6TgREJtouAAgQpiUVJnBIwO1nT5BtqDAYSIAscoumSMtmROoM266U76aZzhMBMQuGi4AiCBm8CgpuZq/OJQj+9VUt4EAxBx78IB/oYxSmWcrtcwvewCEFg0XAEQQk15D5vwLArX97APZ4mKHiQDEEmutvH89Lx084B9ISpbvujtkjHEbDIhhNFwAEGHMsEulhCOrhO3dLTvvC7eBAMQMO/tzaXHwUkIz+gaZRk0dJgJiHw0XAEQYU7ueTN8BgdpO+a+s5zlMBJw6a63s4UOy2Ttld2yR3fK97OYN/vpwrmzp/k8IO7tnl+x/XgkOtO/iv4wZQFix0QIARCBzwRWy33zu34x0+2Ypa57Uta/rWMBx2bzD0obVsmu/ld20TsreKe3dLeUdPvaLjE+qXVdq0ESmYROpcTOZ1u38//fxe+FQsSUl8v7+tJTvX5VQqWny3fQr5hioAjRcABCBTKOmUre+0qLZkiRv8vvydenDfRaIOHbXNtmFs2UXfSNtXOv/JUGl3sCT9uyS9uySXbnYPyRJqWlS63YyZ3eV6dJHpl7DkGePJ/ajt6Q1ywO1ueYWmdr1HCYC4gcNFwBEKN+FV8o70nBpwxr/D0ttO7kNBUiy+Xmys2fIfvmptHnDyb8wIVFKTPIfH9n/6ZjycqXlC2WXL5R993X/Ga+ufWX6DpTq0ShUhl0yV/aT9wO16XWeTJ+B7gIBcYaGCwAilGl5ltS+i7RqqST/Wa4EGi44ZPfslv18guyXU/0NUUWqp0mt2sm0aus/U1ungVS3vpSWUe7yNVtSIh3OlXJzpN07ZHdslXZskf1+rb+Jsz+4b3HbJtltm2Qnj9Oes9orZcBFsh17yqSlh/Ezjn5213b/pYSlGp0pc8MvOVsOVCEaLgCIYL6LrpR3pOHSysWyG9fJNG/tNhTijj2YIzt5nOzMyVJF2xQ0aS7Tvb9M975S4+YndV+QSUiQMmr4/2vYVKZTz+DHyz8srV8j++1S2aXzpW2byr226LtVKvpulZScLNN7gMygETLN+HPxQ7agQN6LjwWb42qp8t12n0xKqttgQJyh4QKASNaus9S8jf/eGEn2k/dlbv2t41CIF7aoUPbTD2Q/HR9cbKFUSqrMOUNlBlwo0+jMkH5ck1Jd6tBVpkNXafSNsru2yy6ZIztnZvlLGAsLZb+a6t8gvE0H+S66QurUk7M3kqxXIu+1/5O2BOfL3HhnyL9WAE6MhgsAIpgxRr4RV/p/Sy3JLvpGduc2mQaNHSdDrLOrl8t78zlp17byD5xRV+aCy2X6D5FJrV4lWUz9RjLDL5eGXy675XulLJ2jvM8ny+bsDz5p7Up5z66UmraUuegKmZ7nyPgSqiRfpLHWyr7zmrSkzH5bQy+Rr9e5DlMB8YuGCwAiXde+UsMm0o6tkrWyn46XueGXrlMhRtnDh2Tf/6fsl5+Vf6B6usyIq/yX7yVXcxNOkmnaQjW69VLGNbdo15SPZGdM8i8qU2rLBtlXn5D98N8yF14h02+wTFKSs7wu2M8+9M9Lqc69ZK78ibtAQJxj8wUAiHDG55O58IpAbb/5XHbfHoeJEKvs2lXyHrqrfLNlfP6zI4++It8FlztttsoyydXk6zdICfc/Id99f5W69C7/hN07ZP/1vLz7fy5v2keyhQVuglYxb85M2ff/ERxocZZ8t/zGf88cACdouAAgCpg+A6Qz6vqLkmLZaR+5DYSYYj1P3ifvy/vf3/k3Ki7VtKV8v/tf+a6+WaZ65K4GaFq3U8Iv/yDfA3+T6T3Av5lyqf17Zd99Xd79t8R84+XNmSFbdkXCeg3lu/P/yVRLcZYJAA0XAEQFk5gkM/zSQG1nfSqbe9BhIsQKm3tQ3t8ekh3/puQdWYrd55O57Dr5fv9//u0JooRp2kK+n/+PfA+/KHP+hVJimTsnDuwr03h9HHONl/f1dH+zVbqcfnoN+X71gEyNWi5jARANFwBEDXPeBVJahr8oyCt/jwZwCuyOLfIeuUdasTg4WLuefPc+Jt/IMTKJ0Xmrt6nfSL7rb5fv0ddkhl4iJSUHHzywT/bd1+Td/wt50yfEROPlfTFF9p/PSNb6BzJqynfPX2QaNnEbDIAkGi4AiBqmWorM4FGB2k6fIFuQ7zARoplduVjeI7+Rdm0PDnbrK98fn5Fp3c5dsBAytWrLd/XN8j3ySgWN117Zd14NNl5Fhe6CniLrlcj77xuy/3ohOFjzDPl+84hMk+buggEoh4YLAKKIGTxSKr0f49BB//5DQCV5X0yR98xDwQ1xJf8lhLf9TiYtcu/VOlXlGq8hF0uJZVYtDDRet8ibPjFqGi+bd1je84/ITvlvcLBWHfnueYS9toAIQ8MFAFHEpNfwX1p4hP3sA9niIoeJEE2stfImjfOfESm9Xyu5mny33ee/hDDGNww2tWrL96Ofy/doBY3X/r2y77zib7w+j+zGy27fLO+xe6Ws+cHBJs3lu/dRLiMEIhANFwBEGTPsUinhyL01e7Nl533hNhCigrVW9r2/y3747+BgrTry/fYxme793QVzwNSqc/zG6z9HGq9Px8sezj32G1Ux65XI+/QDeX+6W9q2KfhAl97y3fe4TL2GzrIBODYaLgCIMqZ2XZm+AwO1nTJetvRsBVAB65XIvvE32allthOo39j/Q3qz1u6CORZovB55xX9/5A8br/f/Ke/en8p79zXZ3TvcBZVkt2+R99ff+ffYKnNW21x0pXy33y+TUt1hOgDHE53LDwFAnDMXjpb9Zrp/VbLtm6Wl86RufV3HQgSyXonsP/8mO3tGcLBZK/nuepAlw48wZ9SR+fEtshdeIfvJ+7JffioVF/sfLMiTnfax7LSPpQ7d5BtwgdS5d5Wt4Gj3ZstOfEf262nBy0Al/+IY198h88MNnwFEHBouAIhCpmFTqVs/adE3kiTvk/fl69on5u/BQeVYz5N947nyzVbm2fLd8QeZ6mnugkUoc0YdmWt+ITviStnPJ8nOmiIdPhR8wsrF8lYuljJqyvQ8V6b3eVKrdjK+0F8wZLdvkZ31iT/DD+7TNL0HyFxzi0zpNhEAIhoNFwBEKd+IK+Udabi0YY20epnUrrPbUIgY1vNk33zOfya0VNtO8t35/2RKV7pEhUytOjKjb5AdcZXsN9P9e97t2Bp8wsEDsjMm+cdr15Xp2EOmQzepXefTWuXR5uyTzVrgP5u1dtXRT6jXUL6rfirD2WwgqtBwAUCUMs3bSB26SiuXSPKf5Uqg4YKOLJDxziv+H9xLZZ5Ns1VJJiVVZvAo2UEjpTXLZWdNkV00WyopDj5pb7bsF5/KfvGpZHxSk2YyZ7aUmraUadxMqnmGVLOWlF7D/7jnSV6JdGCftHuHbPZOafN62W+X+S8Prkit2jKjfiRzztCo3YwaiGf8qQWAKOa78Ap5RxourVwiu3GtvxFDXLMT/iM7Y3JwoE0H+e78I83WKTLGSG07ybTtJJt7UHbRbNn5X0rfLpNsmfuqrCdt+V52y/eSZsie7gdu2lLm/OEy/YfKVKt2uu8GwBEaLgCIZu06Sy0z/ZcU6shZrlvvcxwKLnnTJ8pOeCc40DJTvrv+KJOS6i5UDDFpGTLnDZfOG+6/BHD5Iv8vO1YukQ4eOP0PULuuTOdeMucOk5q15r5MIAbQcAFAFDPG+M9yvfiof2DRbNkdW/yLaiDueHNnyb7zSnCg0Zn+M1ssGR4WpsYZMv2HSP2H+Ldm2L5ZdvN6afMG2U3rpeydUs4+qfAYmyjXqi3VbeDfP6tNB5l2naR6jWiygBhDwwUA0a5rH6lhU2nHFsla2U8/kLnxTtepUMXst1my/3gmOFC7nnx3PySTUcNdqDhifD6pSXOZJs2lvoMC49ZaqSBPOpgjGSP5EiSfT6qeJpPMZYJAPGDjYwCIcsbnk7nwikBtZ8+Q3ZvtMBGqmt22Sd4LjwYXc0ivId+vH5KpXddtMMgYI5NSXaZeQ5m6Dfwbl9eqTbMFxBEaLgCIAabP+VLpD9clxbLTPnIbCFXGHtgn729/kvJy/QNJyfL98g9cVgoAEYKGCwBigElMkhl2WaC2X3wqeyjHXSBUCVtQIO+5h6U9u/wDxsh381iZ1u3cBgMABNBwAUCMMOcNl9Iz/EVBfvllwRFzrLWyb/xN+v67wJi56qcy3fs7TAUA+CEaLgCIEaZaisyQiwO1/XyCbEG+w0QIJztlvH8vqCPMoBEyQy9xmAgAUBEaLgCIIWbQSKnakf2WDh2U/fIzt4EQFnbZQtkP3gwOtO8ic/XPWU4cACIQDRcAxBCTliEz4IJAbT/7ULa4yGEihJrdsVXeq09I1voH6jaQ75bfyCQkuA0GAKgQDRcAxBgz9FIp4cg2i/uyZb/53G0ghIzNOyzv+b8EVySsliLfHb+XSWevLQCIVDRcABBjzBl1ZPoPDtR20jjZIs5yRTvrefJef9K/wfURvp/cLdO0hbtQAIATouECgBhkRo4JnuXau1v266luA+G02Y/flpbOC9Rm5BiZHqxICACRjoYLAGKQqVNf5rxhgdpOek+2qNBhIpwOu/Br2UnjggNdestcco27QACAk0bDBQAxylx0lZR45CzX/j2yX3zqNhBOid22Sd4/ngkONDpTvp+NlfHxTzgARAP+tgaAGGVq15U5/8JAbT95X7agwGEiVJYtyJf30uNS6X5q1dP8i2SkVncbDABw0mi4ACCGmYuulJKS/cWBfbKzJrsNhJNmrZX994vS9s2BMd9Px8o0aOwwFQCgsmi4ACCGmVq1ZQZeFKjtlPGy+XnuAuGk2a+mys6ZEajNhVfIdOnlMBEA4FTQcAFAjDMXXiElV/MXBw/IzpjkNhBOyG5aL/v2y8GBszrIXHadu0AAgFNGwwUAMc7UqCUzeFSgtp9+IJt32GEiHI/NOyzv5cel4iN7p2XUlO/nv5FJSHAbDABwSmi4ACAOmOGXS9VS/UXuQdnpH7sNhApZa2XfeFbatd0/YIx8N4+VOaOO22AAgFNGwwUAccBk1JAZenGgtp99JHv4kMNEqIidMUl24deB2oy6WqZDN4eJAACni4YLAOKEGXaZVLqceF6u7FTOckUSu+E72XF/Dw607yIz6mp3gQAAIUHDBQBxwqSlywy9NFDbaR/JHsxxmAilbO4h/31bJcX+gZq1/ZcS+rhvCwCiHQ0XAMQRM/QSqXq6v8jPk508zm0gyFor7x9PS3t2+QeMT75b7pGpcYbTXACA0KDhAoA4YqqnyYy4MlDbGZNld+9wmAj2sw+lpfMCtbn8OpnMju4CAQBCioYLAOKMGTxKql3XX5QUy374lttAccyuXSk7/o3gQKeeMheMdhcIABByNFwAEGdMUrLMpdcGajtvluzGdQ4TxSd78IC8l/9X8jz/QO268v30bhkf/zQDQCzhb3UAiEOm70CpSfNA7ZU9y4Kws54n77Unpf17/AMJCfLdcq9Meg23wQAAIUfDBQBxyPgS5LvixuDAyiWyKxa7CxRn7OT3pJXB+TZX3iTTup3DRACAcEl0HSAaFRYW6sMPP9TXX3+t7Oxspaenq0uXLrr66qtVp04d1/EA4OR07CG17SStXiZJ8t77u+x5g2US+KchnOyqpbIf/yc40K2vzJBL3AUCAIQVZ7gqqbCwUH/+85/1/vvvKz8/Xz179lSdOnU0c+ZM/fa3v9WOHaz2BSA6GGPku+onkjH+ga0blcdmyGFl9++V99r/SfbIfVv1Gsp3069kSr8GAICYQ8NVSR988IFWr16tzMxMPfPMM/r1r3+tRx55RDfccINycnL04osvuo4IACfNNG8j029woD749qvyDh10mCh22ZISea8+IeXs9w8kJsn3i9/KlO6LBgCISTRclVBcXKwpU6ZIkn72s58pJSUl8NioUaPUvHlzrVq1SuvXr3cVEQAqzVx+nVTN//eZzdmvQ+/9w3Gi2GQ//o+0ZnmgNj/6uUzz1g4TAQCqAg1XJXz77bfKzc1VgwYN1LJly6Me79OnjyRpwYIFVR0NAE6ZqVVH5qLgZsiHJ70nu3Obw0Sxp2DhbNnJ4wK16T1A5vwLHCYCAFQVGq5K2LhxoyRV2GxJUqtWrco9DwCihRl2qVS7nr8oLpb33t/dBoohJbt3av/TDwUHGjaVuf527tsCgDhBw1UJ2dnZknTMlQhr165d7nkAEC1McjWZK28KDiydJ5s131meWGGLi7T/iT/IHjzgH0hOlu/W38qkpLoNBgCoMqz9Wwn5+fmSpGrVqlX4eOk9XaXPO5GxY8ceNZacnKzHHntMklS3bt1TiRkSiYn+b4169eo5yxBrmNPQY05Dy150ufZ9M02Fy/37Q5lxr6vuuYNlqqWc4JU4lpxX/k+HVwfv26p5671K7drTYaLYwJ/90GNOQ485Db1onVPOcFWCtfa0HgeASGaM0Rm33iv5EiRJJTu36dB/33ScKnrlzZyiw5PfD9SpQ0YpdfBIh4kAAC5whqsSUlP9l4AUFBRU+HjpeNnVC4/nySefPO7j2dnZzpq40t8c7N6928nHj0XMaegxp6FXr2kLVb/kah3+8G1JUu74fymvcx+ZBo0dJ4sudtN6eS88GqgTW7VVwegb+V4NEf7shx5zGnrMaei5mlNjjBo1anTKr+cMVyWUXuK3Z8+eCh/fu3dvuecBQDRKv/pnUq0j96oWF8v7z8ucwa8Em3tQ3ouPSoWFkiSTUUO1fvuoTHLFl6MDAGIbDVclNG/eXJK0YcOGCh8v3X+r9HkAEI18qdXlu/pnwYEVi2Xnf+kuUBSxnifvtSel7J3+AWNUa+yflNjg1H8zCgCIbjRcldCuXTtVr15dO3furLDpmjt3riSpe/fuVR0NAEKrxzlSh26B0v7nleBKezgmO+EdafnCQG0uvVbVuvVxmAgA4BoNVyUkJibqwgsvlCT9/e9/L7ca4cSJE7Vx40a1a9dObdq0cRURAELCGCPfdbdJpZfBHcqRfec1t6EinF06X3biO8GBrn3KbSgNAIhPLJpRSaNHj9ayZcu0evVq3XXXXWrXrp2ys7P13XffKSMjQ7fffrvriAAQEqZeQ5nLr5N993VJkp03S7b3+TJdejlOFnnsrm3yXi+zEFKDJvL95G4ZH7/XBIB4x78ElZScnKwHHnhAV1xxhZKTkzV//nzt2rVLAwYM0OOPP66GDRu6jggAIWMGj5JatwvU3r9fkD2c6zBR5LH5h/0rEuYdmZdqKfLd9juZ6mlugwEAIgJnuE5BcnKyrr76al199dWuowBAWBlfgnw33invT3dJxcXS/j2y416Tueku19EigvVK/ItkbN0YGDM3/kqmSTOHqQAAkYQzXACA4zKNzpQZ9aNAbb+eLrvwG4eJIocd/y9p6bxAbS4YLV+vcx0mAgBEGhouAMAJmQtGSy0zA7X35nOy+yrekzBeeN9Ml/10fHCgS2+Z0de7CwQAiEg0XACAEzKJifLdPFaqluIfOHxI3j+elvU8t8EcsWtWyP7r+eBAk+by3TxWxpfgLhQAICLRcAEAToqp31jmRz8PDqxaKjvtI3eBHLHbN8t7/mH/PW2SlFFTvjv/n0xKdbfBAAARiYYLAHDSzDlDpe79A7Ud/y/Zdd86TFS17P698p55SCpdqTExSb7b75epU99tMABAxKLhAgCcNGOMfDfcIdWq4x8oKZb30uOyOfud5qoKNv+wvGf/LO3Z5R8wRr6b/0emTXu3wQAAEY2GCwBQKSYtQ75bfiMlHLlfaf8eea/8r2xJidtgYWSLCv17bW1aFxgzY34m06P/cV4FAAANFwDgFJizOshc9dPgwOplsh+86S5QGNmSEnmvPCGtWhoYM0MvlW/oJQ5TAQCiBQ0XAOCUmMGjZHoPCNT20w/kzZ3lMFHoWc+T/ecz0pI5gTHT+3yZq37iMBUAIJrQcAEATokxRuaGO6QmzQNj9p/PyK5e7jBV6FhrZd9+SXbOzOBgl94yP7lbxsc/nwCAk8O/GACAU2aqpch3+++k9Az/QHGxvBf+Irttk9tgp8l6nuy/X5CdNSU42K6zfL+4VyYx0V0wAEDUoeECAJwWU7+xfL/8f1JSsn/gcK68Zx6S3b/XbbBTZL0S2X/+TfaLT4ODLTPlu+P3MqWfIwAAJ4mGCwBw2kzrdvLdPFYyxj+wd7e8px+QPXjAbbBKsiUlsq8/JTv78+Bg63by3f2QTEqqu2AAgKhFwwUACAnTvb/MmJ8FB7ZulPd/f4iapsvmH5b33J9l530RHMzsKN/dD8pUT3MXDAAQ1Wi4AAAh4xt6icyIq4IDpU1XhG+MbPftkffX30nLFwUHO3SV71cPyKRUdxcMABD1aLgAACFlLrtOZsSY4MDWjfKe+L3snt3uQh2H3bJB3qO/kTZvCA527y/fL/8gU62au2AAgJhAwwUACCljjMxl18qMLNN0bd8s75H/kd2wxl2wCnjfTPc3W/uyA2Nm+GX+1QhZIAMAEAI0XACAkDPGyFx6rczFPwoO5uyX97/3y5v/lbtgR9jCAnlvPCv7j2ekwkL/oPHJXPML+a76KftsAQBChs1EAABhYYyRueQaeWfUlX3rRamkRCoqlH3lr/K+X+O/9NDBWSS7cZ28fz4jbfk+OJiWId/Pxsp06lHleQAAsY2GCwAQVr7zhsvWbyTvhUelw4ckSfazD2VXLJbvp7+WadaqSnLYggLZCW/LTv1I8rzgAy0z5fvFb2Xq1KuSHACA+MI1EwCAsDNtO8l3/xNSw6bBwa0b5T1yj7yP35bNzwvbx7aeJ7vwG3kP3Sn76Qflmi0z5GL57n2UZgsAEDac4QIAVAnToLF8f3hK9oM3ZadP8A+WFMtOeEd21hSZUVfLnDdcJjEpJB/PWitlLZD30b/Lr0AoSfUaynf9HTLtu4TkYwEAcCw0XACAKmOqVZP50c9lO/eS98+/BVcHzNkv+/bLsp9+IHPuMJn+g2Vqn9pZJ5uzX3buLNmvp0lbN5Z/0OeTGXaZzMU/Zsl3AECVoOECAFQ506GrfA8+Kzvlv7LTPw6uFLhnl+xHb8l+/LbUvotMh24yrdtKzdscc4ENW1QkbVwru3al7Orl0qol/gU6fqhbX/kuvVamSfOwfV4AAPwQDRcAwAlTPU1m9A2yg0f6Lyv8amrw/iprpZVLZFcukZWkhESpdl2perpUPU3y+aSDOdKhA9KB/VJJ8bE/UKee/kareesq+KwAACiPhgsA4JSpVUfm+jtkR46RnT1D9pvp0q7t5Z9UUizt3nHyb3pGXZl+g2T6DZZp2CS0gQEAqAQaLgBARDC168mMHCM74irpu5WyyxbIrl8tff+dVFhw/BcnJ0stMmXadJBp10lq21HGl1A1wQEAOA4aLgBARDHGSJlny2SeLUmyxcXSjs3+hTVyc/17eZUUSxk1ZTJqShm1pAaNQra6IQAAoUTDBQCIaCYxUWra0n/sOAsAAJXFxscAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJomuA+DYjDGuI0REhljDnIYecxp6zGnoMaehx5yGHnMaesxp6FX1nJ7uxzPWWhuiLAAAAACAMrikEAAAAADChIYLFbrvvvt03333uY4RU5jT0GNOQ485DT3mNPSY09BjTkOPOQ29aJ1T7uFChQoLC11HiDnMaegxp6HHnIYecxp6zGnoMaehx5yGXrTOKWe4AAAAACBMaLgAAAAAIExouAAAAAAgTGi4AAAAACBM2IcLAAAAAMKEM1wAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECY0XAAAAAAQJjRcAAAAABAmNFwAAAAAECaJrgOgauTn52vevHlau3atvvvuO23cuFHFxcW65pprdNlllx33tXv27NG7776rpUuX6tChQ6pbt6769++vyy+/XMnJyZXOsmXLFo0bN04rVqxQfn6+GjZsqEGDBmnEiBHy+aL/dwAPPvigVq5cedznGGP07rvvntT7zZw5Uy+88MIxH+/fv7/uvvvuykSMOitWrNBDDz10zMfPOuss/eUvf6n0+y5cuFAff/yxvv/+e0lSixYtdMkll6hHjx6nGjVqbN26VfPnz1dWVpa2b9+uAwcOKC0tTW3bttXIkSPVvn37Sr1fPH2fFhYW6sMPP9TXX3+t7Oxspaenq0uXLrr66qtVp06dSr1Xbm6u3nvvPc2bN0/79+9XrVq11KtXL40ZM0ZpaWlh+gwiR0FBgZYuXaqFCxdq3bp12r17tzzPU8OGDdWnTx+NGjVKKSkpJ/1+d9xxh3bv3n3Mx5966ik1adIkFNEj2on+Hbr//vvVtWvXk36/eP8+PdG/QaXGjBmjK6+88oTPi6fv0/Xr1ysrKyvw8+e+ffuUlJSkt95667ivmzVrlqZMmaItW7YoMTFRmZmZGj16tNq2bVvpDJ7n6ZNPPtHnn3+uHTt2KCUlRWeffbbGjBmjpk2bnuqnVik0XHFix44deu65507pdX/4wx+Uk5OjM888U+3atdP69ev13//+V8uWLdMDDzygpKSkk36/NWvW6M9//rMKCgrUpk0b1atXT6tWrdKbb76p1atXa+zYsTLGVDpnJOnatavq1atX4WPr16/X5s2b1a5du0q/b/PmzdWiRYujxs8666xKv1e0atCgQYVz16BBg0q/1+TJk/XPf/5TCQkJ6tSpkxITE5WVlaXHH39cN910k0aMGBGKyBHrz3/+s/bu3avU1FSdddZZyszM1JYtWzRv3jzNnz9fN9xwg0aOHFnp943179PCwkL9+c9/1urVq3XGGWeoZ8+e2r17t2bOnKlFixbp4YcfVsOGDU/qvQ4ePKg//OEP2r59uxo0aKBevXppy5Yt+uSTT7R48WL95S9/UUZGRpg/I7e++uorvfzyy5KkM888U126dFFeXp7WrFmjcePG6euvv9aDDz6omjVrVup9BwwYUOF49erVTztzNOnTp0+FDWvt2rVP+j34PpVq1ap1zO8pz/P05ZdfSlKl/22Ph+/T999/XwsWLKjUa9544w1NmjRJycnJ6ty5s4qKipSVlaWlS5dq7Nix6t2790m/l7VWTz/9tObMmaO0tDR1795dBw8e1Ny5c7Vo0SI98MADVfLvEw1XnEhJSdHgwYPVpk0btW7dWnPnztX48eNP+LoXX3xROTk5uuiii/STn/xEklRSUqKnnnpK8+bN0wcffKAxY8acVIaSkhI9++yzKigo0A033KBRo0ZJ8p99e/jhhzX3/7d370FRlW8cwL/LgiIIwoAgCCliiISCoZiFpaaOF9LC21RSpk45Xsbt3qSNWf5Sc3KylJo0UvEy2cVSIxqaaJIKMRCDRG4CgeAiBEILyLLs7w/mnFzYZXdhD7DL9/OXnMvLu8eHd9/nnPe874UL+PnnnzFz5szuf9B+oKsnhq+//joA4MEHHzS7XOFu4kAWHByMDRs29LiciooKJCQkwMHBAdu2bUNQUJC4/Y033kBCQgImTZoEHx+fHv+u/srPzw8rV67EfffdB3v7/74KkpOTcfDgQSQkJCAsLMzsu3+2HqenT59GXl4egoKCsHXrVrEze+7cORw9ehQfffSRSXfCgfZORWVlJSIjI/H8889DLpcDAOLj45GUlIQjR45g48aNkn2W/sDe3h5z587FwoULdf7eamtrsWvXLhQXF+Pw4cPYvHmzWeVaop2wBbGxsfDy8upRGYxTYOTIkQZj6tKlSzh//jw8PDwQEhJiVrkDIU6DgoIwevRoBAYGIjAwEM8++2yXx+fk5OC7776Di4sLduzYIbYL+fn5ePPNNxEXF4eQkBAMHTrUpN+fkpKCtLQ0+Pj4YPv27XBzcwMApKWlYe/evfjggw/w/vvvi3EtFesfv0UmGTFiBNatW4fZs2cjICDApKF7hYWFyM3NxbBhw7By5Upxu1wux9q1ayGXy/H999+jtbXVpDqkp6dDqVRi1KhRYrIFtCeDa9asAdDeabFVlZWVKCwshIODA6ZNm9bX1RnQEhMTodFoMGfOHDHZAgBfX1889thj0Gg0+P777/uwhtLbunUroqKidJItAJgzZw7CwsLQ1taG33//vY9q1z+1trYiKSkJALBmzRqdJwfR0dEYNWoUcnNzce3aNaNl1dXV4fz58zrtqSA2Nhaurq5ITU1FXV2dxT9Hf/LQQw9h7dq1nW5uuLu7i98L6enpJn/PkGUxTo0Tnm5Nnz7dJl6LsLRHH30Uy5cvR0REhJjsdOXs2bMAgJiYGJ12ISgoCHPmzEFjYyNSUlJM/v1Cv/LJJ5/U+f333XcfJk+eDKVSiYsXL5pcXncxMsigzMxMAEBERESnYYNubm4YP348VCoV8vLyTCovIyMDQHuQdxQQEABvb2+UlZWhqqqqhzXvn4RGOSIiwqaGC1gjIbb1xaKQDAvxOhCNGjUKQPtTBvrP1atXoVKp4O3tjYCAgE77p06dCgAmDZ+5dOkStFotQkJCOnVCHBwcEBERgba2NmRlZVmi6lZJiEO1Wo2GhoY+rs3AxDjtWnNzs9hZnz59eh/Xxvq1tLQgJycHgP7vZ2Gbqd/PVVVVKC8vx6BBg3Dvvff2uLye4JBCMkiYSEBfx0LYnpOTg9LSUtxzzz1GyystLTVanlKpRGlpaY+HQPRHQsLVneGEQPv7XwkJCWhqaoKbmxtCQ0PNHr5g7W7cuIETJ06goaEBLi4uCA4ORnh4uFl3FVUqFaqrqwFA77tGHh4ecHFxwc2bN9HY2Dggk2OlUgkAJt2N7MiW49RYGzZmzBid43pSVkBAAFJSUsR2eCAS4lAul5s8fEhw5swZ3LhxAw4ODvD390dkZCRcXV2lqGa/9tNPP+Hff/+FTCaDj48PIiMj4enpafL5jNOupaen4/bt2wgICIC/v7/Z5zNOdVVUVECtVsPV1VXvBERCHJrSxgL/9WP9/f07jeboTnk9wYSLDKqpqQEAg7NuCduFzqsxwnGGyhNe4jW1PGuSn58PpVIJFxcXs2aGulNmZqb4ZAZofxE1JCQECoWiWx1ja5SXl9fpiepdd92FF1980eT3rYT4cnZ2Njj7mYeHBxoaGlBdXY277rqrZ5W2Mjdu3BDjbPLkyWafb8txask2TDjG0OQF5ravtigxMRFA+0RE5kzOBADHjh3T+fnIkSN45plnMGvWLIvVzxp0fFc7ISEBS5YsMWkmPYBxasydwwm7g3Gqy1gb6+joCGdnZ6hUKjQ1NWHIkCE9Kq8345cJFxnU3NwMAAanfh88eLDOcaaWJ5zXkdD5NbU8a/LLL78AaJ8aW99dlq64ublh2bJlmDJlCry8vNDS0oLCwkIcP34cV65cwa5du/DOO+/Y9NhxJycnLFq0CFOnThUTq5KSEpw8eRIFBQXYsWMH9uzZY9LTKGNxeOc+W4zFrmg0GsTFxUGtVuP+++8Xn9iYYiDEqSXbMGNlCdtv375tdj1tQWZmJlJSUiCXy7FixQqTz4uIiEBoaCjGjBkDV1dXKJVKpKSkIDExER9//DGGDh1q1gxn1mr8+PGYNWsWxo0bB3d3d1RXVyMtLQ1ff/01Tp06BScnJ5NmYmWcGlZXV4fs7GzY2dkhKirKrHMZp/oZ63cC7TGnUqnQ3NxsNOGydD+2J5hwWYn33nsPZWVlZp2zceNGjB07ttu/U6vVAoDBadqF/ZZi6fK6y9LXurW1VZx8oDvDCcPDw3Weijk5OWHy5MkIDQ3Fq6++imvXruG3334zu8HvTT29pgEBAZ2GtISGhuLtt9/G9u3bkZubi6SkJMTExBgt11hcWwsp2oT4+HhcvXoV3t7eWLt2rVll20KcGmOsjTKnDbOVOJRCeXk5PvzwQ2i1WsTGxuod+mvI6tWrdX729/fHU089BV9fX3zyySc4fvz4gOjIdkxSfX19ERMTg8DAQPzvf//DqVOnMHv2bKNraTJODUtNTUVbWxvCw8PNfnrPONXPlHjrTl+xP8QvEy4rcfPmTVRUVJh1Tk/vOAl3DgyV09LSAgAmL0rp6OgIlUplsDxhuzmLXErB0tc6KysLDQ0N8PHxsehaD46Ojpg/fz7i4+ORlZXVrzuyUsWvnZ0dFi9ejNzcXFy+fNmkhEuI667uaPWXWOyKpa/pl19+ieTkZAwbNgxbtmwx+50ZQ6wpTo0x1iaaEzfG4lAoq6snsbaopqYG77zzDlQqFaKjoy22Ht6sWbPw+eefo7KyElVVVTb5nrApwsLCEBgYiKKiIuTn5yM0NLTL4xmnhvX0vWx9BnqcGmtjAfP6nsIx/aHfyYTLSuzatavXf6eHhweKi4vFd7k6Erab+gKup6cnVCoVampqxNmn7vTPP/+YVZ5ULH2theGEUsxgJAyv6+9T8koZv8Iis6ZeAyG+hCEJ+hpac2O7L1jymiYlJYnDjLZs2WLywr2mspY4NUaIB0NtojltmHCMcE5H1hCDllZfX48dO3aguroaM2bMQGxsrMXKtrOzg7e3N27duoXa2toB15G904gRI1BUVGTS3yPjVL/y8nIUFxfD0dERU6ZMsVi5Az1OjbWxzc3NUKlUcHZ2Njqc0JTyejN+rXcwPUlOGMZRXFysd7+w3dRJBYQky1h5+pIxa9XY2ChONypFwqVSqQD07ycxUjP3Gjg7O4uNq76ZtWpqatDQ0ABPT88BMUPh+fPn8dlnn2Hw4MF47bXXzBq+ZSpbiVNjbZiw/pYpbdhAbA+70tTUhJ07d+L69euIjIzEunXrLD4MyFbisKfMuQ6MU/2EG6mRkZEWf7o3kOPU19cXDg4OqK+v15skmdvvFL7PysrK9K7lZ255PcGEiwwS1izIyMiAWq3W2VdXV4fc3Fw4OTkhODjYrPLS0tI67SsuLoZSqYSfn59N3dFJS0uDWq3GuHHj4O3tLUn5AMya3MDWXLhwAYDhaYv16SoWhfft9K3ZYWsyMzMRFxcHuVyOl156yeS/ZXPZSpwGBwfDyckJSqVSbwdUiEVTYic8PBwymQy5ubm4deuWzj61Wo2MjAzIZDJMmjTJMpXvx9RqNd59910UFRUhLCwMCoXC4pOrlJWVoaKiAoMHD8bIkSMtWrY1qa+vR25uLgDT2kzGaWdarRa//vorAMsOJwQYp4MGDRKHuer7fha2RUREmFSel5cXRo4ciZaWFp3Zc7tbXk8w4SKDxo4di3HjxuHWrVs4fvy4uF2j0eDQoUPQaDSYN29ep1n39u/fD4VCgfT0dJ3tkZGR8PLyQmlpqbjyN9D+iPjTTz8FAERHR0v4iXqfOWO8FQoFFApFp6EbiYmJncbPt7a24osvvkBaWhoGDRqEGTNmWKzO/VFycnKnhU+1Wi2Sk5Px3XffQSaTYe7cuZ3OM3RNFyxYADs7OyQnJyM/P1/cXllZidOnT8POzs5i7470V1evXsXevXsBtF+nsLAwk84byHFqb2+PefPmAWifYOTOz3vu3DmUlpYiODhYZ2KSpKQkKBQKnDhxQqcsd3d3PPDAA2htbRXbU8GxY8dQX1+PqKgoq59K35i2tjbs27cPf/31F8aPH4+XXnrJ6Eyuhq5pVlaW+JTxTqWlpdi7dy+0Wi1mzZpl9kyx1iY/Px85OTmdJheoqqrCnj17cPv2bUyePFlnqmzGqelyc3Nx8+ZNuLu7d/kOHOO0exYuXAigfUmDyspKcXt+fj5+/PFHDBkypNO0+YWFhVAoFHjrrbc6lSf0K48fP65z0+DChQv4448/4OXlZdFhoYYMzP/NAWrPnj3imG3hUe0PP/wgrpLu5uaGl19+Weec9evXY+vWrUhMTEROTg78/PxQVFQEpVKJu+++W+8kBdXV1aioqEBjY6POdnt7e2zatAlvv/02jh49it9//x2enp64evUqamtrMWXKFKvvkN3pn3/+wZUrV2Bvb49p06YZPV6YAKHjY+/Dhw/jxIkT8PPzg6enJ9RqNUpKSlBbWwsHBwds2rTJ4BoptuKbb75BfHw8/Pz8MHz4cADA33//jaqqKshkMqxatUrv0xND19TX1xcrV67E0aNHsW3bNkycOBFyuRx//vknWlpaxNmibNnu3bvR0tICLy8vXLx4UWwH7hQcHIyHH35YZ9tAj9OYmBhkZ2cjLy8PmzdvRnBwMKqrq1FQUAAXFxesX79e5/j6+npUVFSgtra2U1mrVq1CQUEBLly4AIVCgcDAQJSVlaGsrAze3t54+umne+tj9ZmkpCTx5pyLiwsOHTqk97jY2FhxQVhD1zQ/Px9ffvklhg8fDm9vb7i6uqKqqgrFxcXQaDQICQnBE088Ie0H6gcqKioQFxcHd3d3+Pj4wM3NDTU1Nbh27RrUajX8/f3x3HPP6ZzDODXdnWtvdfUklnHaLjMzE1999ZXOttbWVmzZskX8ecmSJeLIgIkTJ2LBggVITEzEK6+8ggkTJkCj0eDPP/9EW1sbNm3a1GlSp9u3b4uLJnc0c+ZMXLp0Cenp6VAoFJgwYQIaGhpw5coV8bupN5JbJlwDSElJCW7evKmzraamRky+hI7snXx8fLB7926cOnUKWVlZSE9Ph4eHB2JiYhATE2N0StmOxo0bh507d+LUqVO4cuUKSkpK4O3tjejoaCxcuNCq1+jp6Pz589Bqtbj33nt7NOPb0qVLkZ+fj+vXr6O8vBxarRYeHh6YPXs2oqOjbT4xANrvUF2+fBnl5eXIzs6GRqOBu7s7pk+fjvnz53dr+YPo6GiMGDECZ8+eFYfYjBkzBosWLerWgr/WRnhPoKqqClVVVQaP65hwGTJQ4nTQoEHYtm0bTp8+jdTUVFy8eBHOzs546KGHsGLFCrNevnZ1dRXbw4sXLyI9PR3Dhg3DvHnzsHz5covNFNmf/fvvv+K/O46KuNOyZcvEhMuQ8PBw1NTUoKioCKWlpWhsbMSQIUMQHByMqKgozJw506a+YwwZO3Ys5s6di4KCApSXlyMvLw+DBw/G6NGjMW3aNMydO9es727G6X/UarU4DK2772UPtDitr69HQUGBzjatVquzrb6+Xmf/qlWrMHr0aCQlJSE7OxtyuRyhoaFYsmSJ2UPf7ezs8MILLyAxMREpKSnIyMgQJztZsWIF/P39u//hzCDT9pfFj4iIiIiIiGyM7aTQRERERERE/QwTLiIiIiIiIokw4SIiIiIiIpIIEy4iIiIiIiKJMOEiIiIiIiKSCBMuIiIiIiIiiTDhIiIiIiIikggTLiIiIiIiIokw4SIiIiIiIpIIEy4iIiIiIiKJMOEiIiIiIiKSCBMuIiIiIiIiiTDhIiIiIiIikggTLiIiIiIiIokw4SIiIiIiIpIIEy4iIiIiIiKJMOEiIiKygG+++QbLly/H448/jsLCQr3HZGZmYsWKFVi+fDlSU1N7uYZERNQXmHARERFZwOLFizFhwgRoNBrs27cPTU1NOvtra2sRFxcHrVaLBx98EFFRUX1UUyIi6k1MuIiIiCxAJpNh06ZNGDZsGJRKJQ4ePCju02q12L9/P+rr6zFixAisXbu2D2tKRES9iQkXERGRhbi5uWH9+vWQyWRITU3Fzz//DAD49ttvkZ2dDblcjs2bN8PR0bFvK0pERL2GCRcREZEFTZo0CQsXLgQAxMfH45dffsHnn38OAHj88ccRGBjYl9UjIqJeJtNqtdq+rgQREZEtaW1txdatW3Ht2jVxW1hYGF5//XXIZLI+rBkREfU2PuEiIiKyMHt7e6xfv1782cnJCRs2bGCyRUQ0ADHhIiIiksCPP/4o/rupqQklJSV9VxkiIuozTLiIiIgsLCMjA0lJSQCAUaNGQavV4sCBA6irq+vbihERUa9jwkVERGRBwnpbADBjxgxs374dw4cPx61bt3DgwAHw1WkiooGFCRcREZGFtLW1Yf/+/WhoaICPjw9Wr14NJycnbN68GXK5HJcvX8a5c+f6uppERNSLmHARERFZyJkzZ/SutxUUFISlS5cCAE6ePKkzeyEREdk2JlxEREQWUFhYqLPe1pgxY3T2P/bYY7jnnnvQ2tqKffv2obm5uS+qSUREvYwJFxERUQ81NTVh37590Gg0mDhxIh555JFOx9jZ2WHjxo0YOnQoKisrER8f3wc1JSKi3saFj4mIiIiIiCTCJ1xEREREREQSYcJFREREREQkESZcREREREREEmHCRUREREREJBEmXERERERERBJhwkVERERERCQRJlxEREREREQSYcJFREREREQkESZcREREREREEmHCRUREREREJBEmXERERERERBJhwkVERERERCQRJlxEREREREQSYcJFREREREQkESZcREREREREEmHCRUREREREJBEmXERERERERBJhwkVERERERCQRJlxEREREREQS+T8j9AvNQ1hJmwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 960x720 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : x**2 + 10*np.sin(x)\n",
    "fp = lambda x: 2*x + 10*np.cos(x)\n",
    "\n",
    "x = np.arange(-10, 10, 0.1)\n",
    "\n",
    "plt.plot(x, f(x))\n",
    "plt.xlabel('x')\n",
    "plt.ylabel('y')\n",
    "plt.legend([\"$x^2 + 10\\sin(x)$\"])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "757e1afe",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = -2.3787, err = 1.621e+00\n",
      "x1 = -2.4832, err = -1.045e-01\n",
      "x1 = -2.4795, err = 3.666e-03\n",
      "x1 = -2.4795, err = 4.253e-06\n",
      "x1 = -2.4795, err = 5.736e-12\n"
     ]
    }
   ],
   "source": [
    "# 초기값 -4\n",
    "x0 = x = -4\n",
    "\n",
    "Err = []\n",
    "itmax = 5\n",
    "for i in range(itmax):\n",
    "    # Newton-Raphson\n",
    "    xn = x - f(x) / fp(x)\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "cd14fbd5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = -0.2717, err = -1.272e+00\n",
      "x1 = 0.0154, err = 2.872e-01\n",
      "x1 = 0.0000, err = -1.541e-02\n",
      "x1 = 0.0000, err = -2.251e-05\n",
      "x1 = 0.0000, err = -5.067e-11\n"
     ]
    }
   ],
   "source": [
    "# 초기값 1\n",
    "x0 = x = 1\n",
    "\n",
    "Err = []\n",
    "itmax = 5\n",
    "for i in range(itmax):\n",
    "    # Newton-Raphson\n",
    "    xn = x - f(x) / fp(x)\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "b55da5b9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 6.5628, err = 4.063e+00\n",
      "x1 = 4.5472, err = -2.016e+00\n",
      "x1 = 3.0957, err = -1.451e+00\n",
      "x1 = 5.7397, err = 2.644e+00\n",
      "x1 = 4.3537, err = -1.386e+00\n",
      "x1 = 2.5082, err = -1.846e+00\n",
      "x1 = 6.5195, err = 4.011e+00\n",
      "x1 = 4.5492, err = -1.970e+00\n",
      "x1 = 3.1005, err = -1.449e+00\n",
      "x1 = 5.7449, err = 2.644e+00\n",
      "x1 = 4.3563, err = -1.389e+00\n",
      "x1 = 2.5186, err = -1.838e+00\n",
      "x1 = 6.4670, err = 3.948e+00\n",
      "x1 = 4.5496, err = -1.917e+00\n",
      "x1 = 3.1014, err = -1.448e+00\n",
      "x1 = 5.7459, err = 2.645e+00\n",
      "x1 = 4.3568, err = -1.389e+00\n",
      "x1 = 2.5206, err = -1.836e+00\n",
      "x1 = 6.4575, err = 3.937e+00\n",
      "x1 = 4.5495, err = -1.908e+00\n",
      "x1 = 3.1010, err = -1.448e+00\n",
      "x1 = 5.7455, err = 2.644e+00\n",
      "x1 = 4.3566, err = -1.389e+00\n",
      "x1 = 2.5197, err = -1.837e+00\n",
      "x1 = 6.4618, err = 3.942e+00\n",
      "x1 = 4.5495, err = -1.912e+00\n",
      "x1 = 3.1012, err = -1.448e+00\n",
      "x1 = 5.7457, err = 2.645e+00\n",
      "x1 = 4.3567, err = -1.389e+00\n",
      "x1 = 2.5201, err = -1.837e+00\n",
      "x1 = 6.4596, err = 3.939e+00\n",
      "x1 = 4.5495, err = -1.910e+00\n",
      "x1 = 3.1011, err = -1.448e+00\n",
      "x1 = 5.7456, err = 2.645e+00\n",
      "x1 = 4.3566, err = -1.389e+00\n",
      "x1 = 2.5199, err = -1.837e+00\n",
      "x1 = 6.4606, err = 3.941e+00\n",
      "x1 = 4.5495, err = -1.911e+00\n",
      "x1 = 3.1011, err = -1.448e+00\n",
      "x1 = 5.7456, err = 2.645e+00\n",
      "x1 = 4.3567, err = -1.389e+00\n",
      "x1 = 2.5200, err = -1.837e+00\n",
      "x1 = 6.4601, err = 3.940e+00\n",
      "x1 = 4.5495, err = -1.911e+00\n",
      "x1 = 3.1011, err = -1.448e+00\n",
      "x1 = 5.7456, err = 2.645e+00\n",
      "x1 = 4.3566, err = -1.389e+00\n",
      "x1 = 2.5200, err = -1.837e+00\n",
      "x1 = 6.4604, err = 3.940e+00\n",
      "x1 = 4.5495, err = -1.911e+00\n"
     ]
    }
   ],
   "source": [
    "# 초기값 2.5\n",
    "x0 = x = 2.5\n",
    "\n",
    "Err = []\n",
    "itmax = 50\n",
    "for i in range(itmax):\n",
    "    # Newton-Raphson\n",
    "    xn = x - f(x) / fp(x)\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "28a1de15",
   "metadata": {},
   "source": [
    "### 문제점\n",
    "- 수렴하지 않을 수 있다.\n",
    "    * Local extrema 근처에서 진동할 수 있다.\n",
    "\n",
    "- 초기 추정 값이 중요함\n",
    "\n",
    "### DIY\n",
    "Newton-Raphson method 함수를 구성하시오. 아래 판별 기준을 적용하라.\n",
    "- $f(x_i)$ 값이 tolerance 보다 작으면 멈춘다.\n",
    "- 근사 상대 오차 $\\epsilon_a$ 가 tolerance 보다 작으면 멈춘다.\n",
    "- $f'(x_i) = 0$ 인 경우 에러를 발생시킨다 (`raise ValueError`)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "e63424a2",
   "metadata": {},
   "outputs": [],
   "source": [
    "def newton(f, fp, x0, rtol=1e-6):\n",
    "    # DIY\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "418a65e8",
   "metadata": {},
   "source": [
    "## Secant Method\n",
    "Newton-Raphson 방법은 미분 함수 $f'(x)$ 가 필요하다. 미분을 구하기 힘들 경우 수치적으로 근사한다.\n",
    "\n",
    ":::{figure-md} secant\n",
    "<img src=\"https://upload.wikimedia.org/wikipedia/commons/9/92/Secant_method.svg\" alt=\"secant-fig\">\n",
    "\n",
    "Secant Method (from Wikipedia)\n",
    ":::\n",
    "\n",
    "그림과 같이 ${i-1}$, ${i}$ 번째 결과로 부터 기울기를 대신 사용한다.\n",
    "\n",
    "$$\n",
    "f'(x_i) \\approx \\frac{f(x_i) - f(x_{i-1})}{x_{i} - x_{i-1}}\n",
    "$$\n",
    "\n",
    "이 근사 미분값을 사용하는 방법이 Secant method 이다.\n",
    "\n",
    "$$\n",
    "x_{i+1} = x_{i} - \\frac{f(x_i)}{f'(x_i)} = x_i - \\frac{f(x_i)(x_{i} - x_{i-1})}{f(x_i) - f(x_{i-1})}\n",
    "$$\n",
    "\n",
    "### 예제 1\n",
    "다음 함수 $f(x) = e^{-x} - x$ 의 근을 구하시오. Secant Method로 해석해보자.\n",
    "(초기 값은 $x_{-1} = 0, x_0 = 1$ 로 한다.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "3ea29bcc",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 0.6127, err = -3.873e-01\n",
      "x1 = 0.5638, err = -4.886e-02\n",
      "x1 = 0.5672, err = 3.332e-03\n",
      "x1 = 0.5671, err = -2.705e-05\n",
      "x1 = 0.5671, err = -1.620e-08\n"
     ]
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : np.exp(-x) - x\n",
    "\n",
    "xm, x0 = 0.0, 1.0\n",
    "x = x0\n",
    "\n",
    "Err = []\n",
    "itmax = 5\n",
    "for i in range(itmax):\n",
    "    # Scecant Method\n",
    "    xn = x - f(x)*(x - xm) / (f(x) - f(xm))\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    xm = x\n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f9dd57f6",
   "metadata": {},
   "source": [
    "### 예제 2\n",
    "다음 함수 $f(x) = log(x)$ 의 근을 구하시오.  Secant Method로 해석해보자.\n",
    "(초기 값은 $x_{-1} = 0.5, x_0 = 5$ 로 한다.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "40f2c61f",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x1 = 1.8546, err = -3.145e+00\n",
      "x1 = -0.1044, err = -1.959e+00\n",
      "x1 = nan, err = nan\n",
      "x1 = nan, err = nan\n",
      "x1 = nan, err = nan\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_11119/3825093768.py:2: RuntimeWarning: invalid value encountered in log\n",
      "  f = lambda x : np.log(x)\n"
     ]
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : np.log(x)\n",
    "\n",
    "xm, x0 = 0.5, 5\n",
    "x = x0\n",
    "\n",
    "Err = []\n",
    "itmax = 5\n",
    "for i in range(itmax):\n",
    "    # Scecant Method\n",
    "    xn = x - f(x)*(x - xm) / (f(x) - f(xm))\n",
    "    \n",
    "    # Compute difference\n",
    "    Eps = xn - x\n",
    "    \n",
    "    xm = x\n",
    "    x = xn\n",
    "    print(\"x1 = {:.4f}, err = {:.3e}\".format(x, Eps))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "6c3a71e1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      converged: True\n",
       "           flag: 'converged'\n",
       " function_calls: 44\n",
       "     iterations: 42\n",
       "           root: 0.9999999999998863"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "root_scalar(f, bracket=[0.5, 5], method='bisect')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8f611118",
   "metadata": {},
   "source": [
    "### 특징\n",
    "- 미분 함수 없이 계산 가능하다.\n",
    "- 수렴 속도가 Newton-Raphson과 거의 비슷하나 약간 늦다.\n",
    "- 발산할 수 있다."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "eee2ef8f",
   "metadata": {},
   "source": [
    "### DIY\n",
    "Scecant method 함수를 구성하시오. Newton-Raphson 과 동일한 판별 조건을 사용하자."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "ac943e8a",
   "metadata": {},
   "outputs": [],
   "source": [
    "def secant(f, x0, x1, rtol=1e-6):\n",
    "    # DIY\n",
    "    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4e982bb",
   "metadata": {},
   "source": [
    "## SciPy 활용"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "3df0379f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# root_scalar?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "a7dd9b3a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      converged: True\n",
       "           flag: 'converged'\n",
       " function_calls: 8\n",
       "     iterations: 4\n",
       "           root: 0.5671432904097838"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Function and derivatives\n",
    "f = lambda x : np.exp(-x) - x\n",
    "fp = lambda x: -np.exp(-x) - 1\n",
    "\n",
    "root_scalar(f, fprime=fp, x0=1, method='newton')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "86732f69",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      converged: True\n",
       "           flag: 'converged'\n",
       " function_calls: 7\n",
       "     iterations: 6\n",
       "           root: 0.567143290409784"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "root_scalar(f, x0=0, x1=1, method='secant')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}